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ABSTRACT 



The theory, detailed outline of operating and 
algorithm procedures of a continuous time-invariant, 
linear control state variable systems design and 
cinalysis computer program is presented. The program, 
SVS, which is based on Melsa’s LINCON, was modified to 
demonstrate Controllability, Observability, Bode Plot, 
Root locus plot, Nyquist plot, pole placement, 
Luenberger observer design, optimal control design, 
time response plot and some basic matrix manipulations. 
Worked examples with the program output are included. 
Some options give only numeric data output; others give 
both numeric data and high-resolution graphs. The 
software, which is fully interactive, menu driven and 
user friendly is written in Turbo Pascal to be run on 
the IBM-PC microcomputers. All options are presented 
via option menus and the user will be prompted for all 
input parameters. 
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I. LINEAR CONTROL STATE VARIABLE SYSTEMS 

(SVS) 



A. INTRODUCTION TO THE SVS 

SVS Is software for the IBM-PC microcomputer in 
the analysis and design of continuous time, linear 
control systems. These programs are based on the matrix 
mathematics of state variables and were first developed 
by Melsa [Ref. 1] and adapted for batch use at NFS by 
Desjardins [Ref. 2]. 

The original intent of this thesis was simply to 
take Desjardins’ adapted version of Melsa’ s LINCON and 
modify it for the Turbo Pascal computer language. The 
features of this thesis are: 

- User-friendly as possible 

- Menu-driven program. <Q>quit key always returns to 

the SVS main menu. 

The hierarchical menu structure is three levels 
deep at any point. So the user, before selecting 
his/her option, can get help at every menu. 

SVS was tested with several exajnples and is now 
available to any user on the Naval PostGraduate School 
control laboratory PC’s under the SVS directory name. 

B . SYSTEM REQUIREMENTS 

SVS is a large program and requires at least 512KB 
of Memory to run. The prograun will run on any IBM-PC or 
compatible "MS-DOS" computer and requires a standard 
IBM color Graphics Adapter (CGA) or IBM Enhanced 
Graphics Adapter (EGA) card. It will run on either 
monochrome or color monitor, but all menus are 
supported by color, for that reason they are easier to 
work with a color monitor if available. The graphics in 
the program are in high resolution (640 X 200) mode. 
They only appear in the white-on-black. The graphics 
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can be dumped to an EPSON, IBM-Graphics , or compatible 
printer by using the Shift + PrtSc key. 

The program is written in the Turbo Pascal 
Language. Turbo Pascal has a 64K data segment and code 
sizes limitation. Under this restriction, SVS is 
compiled as five executable programs and nineteen 
"Chain" files. The main program name is SVS in the disk 
files. It has extension .COM. The chain files have 
extension .CHN and they are not themselves executable. 



C. FILES ON THE DISK 

A complete list of the files and a brief 
description of these are below. 



( 1 ) SVS.COM { 

(2) INPUT.COM { 

(3) MATRIX.COM { 

(4) PLOT.COM < 
( 5 ) CHANGE . CHN < 

{ 

(6) POLYNOM.CHN { 

< 

(7) CONTROL. CHN { 

{ 

(5) DETERMIN.CHN { 

( 

(9) EIGEN. CHN < 

(10) BODE.CHN < 

( 1 1 )TIMEPL0T .CHN{ 

( 12) 0BSER.CHN { 

< 

( 1 3 ) LUENBERG. CHN{ 

< 

{ 

( 1 4) 0PTIMAL.CHN < 



The executable main menu module > 
The input/chauige menu module > 

The matrix mathematics menu module) 
The graphics menu module > 

Allows user to change input data 
values > 

Calculates characteristic 
polynomial of A matrix > 

Calculates controllability of the 
system ) 

Calculates determinant of the A 
matrix > 

Shows eigenvalues of the A matrix ) 
Gives Bode plot of the system > 
Gives time response of the system) 
Calculates observability of the 
system) 

Design Luenberger observer to 
achieve a closed-loop poles of the 
system ) 

Optimal control design program ) 
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( 1 5) P0LE.C0M 

( 1 6) SAVE.CHN 

( 17) RETRIEVE.CHN 

(18) INVERSE . CRN 

( 19) INPUTDAT.CHN 

(20) HELP1 .CRN 
(21 )RELP2.CRN 

(22) RELP3.CRN 

(23) RELP4.CRN 

(24) NYQUIST.CRN 

(25) RL0CI .CRN 

(26) R00TS.CRN 



( To calculate feedback coefficients 
( to achieve a desired closed-loop 
( poles > 

( Saves data to desired drive ) 

( Gets data from desired drive > 

( To calculate inverse of the A 
{ matrix > 

{ To permit user to enter the input 
( data for whole system > 

( Relp program for SVS main menu 
{ options) 

{ Relp information for the 
{ input/change menu options ) 

( Relp information for the matrix 
{ mathematics menu option ) 

{ Relp information for graphics menu 
( option > 

{ to calculate Nyquist (polar) plot 
( of the system ) 

( Root locus plotting procedure) 

( Shows plant characteristic 
{ equation roots > 



In addition to these files, there are three 
"system" files that are needed to run the program. 
These are 4X6. EON, 8X8. EON and ERROR. MSG files. 

D. STARTING TO THE SVS 

This software package has two diskettes. The first 
step is to get the SVS main menu on your screen. All 
you have to do is follow these steps: 

1 . Turn on the power ( or, if it’s already on, the 
user reboots the computer ) . 



2. Wait for the operating system prompt. It will 
look like one of the following: 

C:_ 

C>_ 

C\>_ 

C:\_ 

C: \>_ 

<C\ >_ 

(or some other letter). 

The prompt may look somewhat different, depending 
on the computer aind how it has been set up. 

3. Type MD SVS and press <ENTER> key. This opens a 
new directory. 

4. Type CD\SVS and press <ENTER> key. This enters 
the new directory. 

5. Insert disk #1 into disk driver A: 

Type Copy A:**.** and press <ENTER> key. 

6. Repeat step 5 for disk § 2 . 

7. Type SVS and press <ENTER> key. 

Now the user is in the SVS main menu and ready to 
work. Make a choice for further step. 
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II. INPUT /CHANGE MENU 



A. INTRODUCTION TO THE INPUT/ CHANGE MENU 

The first step, of course, is to enter the A, B and 
C matrices into the computer as a common input. This 
is the starting point of the program. For our case, 
the general state variable equations are represented by 
the following equations. 

x(t) = A x(t) + B u(t) (2-1) 

y(t) = C x(t) + D u(t) (2-2) 

where A, B, C are given matrices and D is assumed to 
have a zero elements for this program. 



SVS MAIN MENU *»***» 

(I) Input/Change Plant Matrices Menu 

(G) Graphics Menu 

(C) Controllability 

(O) Observability 

(L) Luenberger Observer Design 

(D) Design of Optimal Control 

(P) Pole Placement 

(M) Matrix Mathematics Menu 

(H) Help 

(0) Quit the Program 
Make Your Selection_ 

Naval Postgraduate School Ismail UNLU 



Figure 2.1 The SVS main menu. 



From the opening menu of SVS, shown in Figure 2.1, 
we choose the ” Input /Change Plant Matrices Menu” 
option to bring us to the Input /Change menu shown in 
Figure 2.2. 



INPUT / CHANGE MENU 

(I) Input Plant Matrices 

(C) Change Current Plant Matrices 

(S) Save Plant Matrices to Disk File 
(L) Load Plant Matrices From Disk File 

(H) Help 

(Q) Quit to SVS Main Menu 
Make Your Selection 



Figure 2.2 The input /change menu. 



This input/change menu allows the user to enter 
common inputs to all prograuns. These common inputs are 
the plant or A matrix, the input or B matrix and the 
output or C matrix. The only restriction is the 

dimensions of the matrix. It must be no greater than 
10. This means that the maximum matrix size has to be 
10X10. However, due to a user decision, a dimension 
size not to exceed 6 is required. The reasoning behind 
this was due, in part, to the printer. Since the 
output format to the printer is Ell, this naturally 
limits us to 6 numbers per line. 12 format places are 
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normally considered necessary for good accuracy. For 
systems with order greater than 6, every row of the 
matrix is continued on the second line. After 
attempting this, it was decided the results were 
difficult to read. Otherwise it is appropriate up to a 
1 0 'th order system. 

B. INPUT PLANT MATRICES 

This option is used to initially enter the A, B and 
C matrices of the state variable equations (2-1) and 
(2-2). For this option, select the "Input plant 
matrices" from the input/change menu. The screen will 
prompt the degree of the plant, which is the dimension 
of the A matrix. The maximum acceptable degree is 10. 
Then it asks for elements of the A matrix and so forth. 

Matrices are entered one element at a time 

beginning with 1,1 and continuing across the row of the 
matrix. The next row is then entered, and the process 
continues until all elements have been entered. After 
the matrix is entered, the complete matrix is 

automatically brought to the screen for review and 
possible element changes. If a change to the matrix is 
desired the user simply enters the row auid column 
number to change after every prompt. Then the program 
asks for corresponding matrix elements to change. 
After being prompted, the chainge is entered. A review 
of the matrix is again brought to the screen. The user 
is again prompted for any more possible changes. This 
procedure continues until all changes have been done. 
The sajne A matrix procedures are repeated for the B and 
C matrices. This input data can be used for the 
options without saving to disk file until quitting the 
program . 
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C. CHANGE CURRENT PLANT MATRICES 

This option allows the user to change input data A, 
B and C matrices quickly and easily. The user can 
change the order and elements of the matrices that were 
previously entered. This powerful combination 
facilitates both input correction and changes, 
especially for higher order matrices during the design 
process. First Figure 2.3 appears on the screen. The 
user can choose one matrix at a time for correction. 
Then the program gives the order and elements of the 
selected matrix. In the beginning the user can enter a 



****** Change Current Plant Matrices Procedure 



Which matrix do you want to change ? PLANT (A) A 

INPUT (B) 
OUTPUT (C) 



Press <ESC> to change it!, 

Then input your choice with <ENTER> key 



Figure 2.3 Change current plant matrices selection. 



correction to the order of the matrix. Then, under the 
new dimension, can make corrections on elements of the 
matrix. The program shows the corrected results on the 
screen. At the end of the program, the user 
automatically returns to the input/change menu. If the 
user wants to change more than one matrix, he/she must 
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choose the "Change Current Plant Matrices" option two 
or more times. A basic example for this option is 
illustrated in Figure 2.4. 

D. SAVE CUERENT MATRICES TO DISK FILE 

This procedure is used to store plant matrices to 
the hard disk or floppy disk file. First it prompts 
drive C as a saving drive. Then asks the user for the 
drive designator (A through E), and filename for the 
problem to be saved. Eight characters of a MS-DOS 
filename are allowed; the program gives a filename 
extension of ".SVS" to each data. This extension is 
used to limit the disk search for appropriate files. A 
drive and filename are supplied by the program, which 
opens the file and stores the data. The procedure 
stores the data as a text file. A text file consists 
of ASCII characters, and is usually designed to hold 
readable Information [Ref. 3]. 

E. LOAD PLANT MATRICES FROM DISK FILE 

The procedure first asks the data drive where the 
problems are stored. After this is done, the program 
calls another procedure called "Directory" [Ref. 4]. 
This procedure uses MS-DOS function, calls and shows 
all available data files on the screen. The user can 
choose one of the files by moving arrow keys with 
<RETURN> key. Then the program opens the file and 
reads it. 

The directory displays only the disk files with the 
extension ".SVS". This eliminates the possibility to 
read other files. 



The A matrix is : 

0.0000E+00 1 .0000E+00 0.0000E+00 0 . 0000E+00 

0.0000E+00 0.0000E+00 1 .0000E+00 0.0000E+00 

0.0000E+00 0 . 0000E+00 0.0000E+00 1.0000E+00 

0.0000E+00 -1.5000E+01 -2.3000E+01 -9.0000E+00 

The order of the system is:4. Change ? (Y/N) 

The order of the system is: 3 

The A Matrix is : 

0.0000E+00 1 .0000E+00 0.0000E+00 

0 . 0000E+00 0.0000E+00 1 .0000E+00 

0.0000E+00 0.0000E+00 0.0000E+00 

Do you want to change any element ? (Y/N) 

Input row to change : 1 

Input column to change : 1 

A( 1 , 1 )= 10 

The A Matrix is : 

1.0000E+01 1.0000E+00 0.0000E+00 

0.0000E+00 0.0000E+00 1.0000E+00 

0 . 0000E+00 0 . 0000E+00 0 . 0000E+00 



Do you want to chajige any element ? (Y/N) 



Figure 2.4 Change current plant matrices example. 
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III. GRAPHICS MENU 



A. INTRODUCTION TO THE GRAPHICS MENU 

From the opening menu of SVS, shown in Figure 2.1 , 
we choose the ” Graphics Menu " option to bring us the 
Graphics menu shown in Figure 3.1. The Graphics menu 
contains five basic options. These are the time 
response plot, Nyquist plot. Root locus plot, 
characteristic equation roots and Bode plot. This menu 
is also supported with the "Load plant matrices from 
disk file" option. This selection allows the user to 
get data from the disk file quickly instead of going 
via the SVS main menu route. 



GRAPHICS MENU “*»»» 

(L) Load Plant Matrices From disk File 

(C) Characteristic Equation Roots 

(B) Bode Plot 

(N) Nyquist Plot 

(T) Time Response Plot 

(R) Root Locus Plot 

(H) Help 

(Q) Quit to SVS Main Menu 
Make Your Selection 



Figure 3.1 The graphics menu. 



B. BODE PLOT 

Bode plot analysis can be accomplished with the 
Graphics menu by selecting the <B> Bode plot option. 
This selection brings Figure 3.2 to the screen. There 
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are two selections for the frequency plot: the 
open-loop Bode plot and closed-loop Bode plot. Input 
data for this option is entered with the Input/change 
plant matrices menu selection which is explained in 
Chapter II. If the user wants the closed-loop Bode 
plot, the program automatically calculates for the 
negative unity-feedback condition. 



****** Bode Plotting Parameters 

Open (0) or Closed (C) Loop Plot? 0 

What is the first frequency to be plotted? .1 
(Example: . 01 , 1 , 100, etc. ) 

How many decades do you want plotted? 4 



Figure 3.2 Bode plot parameters selection. 



The user also must enter the starting frequency and 
number of decades for the plotting. The upper 
frequency limit is calculated based on the number of 
decades. That is, if user selects .1 as the starting 
frequency with 4 decades, then the upper frequency will 
be 1000 rad/sec. 

Bode plot displays two plots at the same time. 
These are plots of magnitude and phase versus radian 
frequency. Magnitude is converted to the decibels unit 
using the relation 

Magnitude^iB = 20 log^ 0 ( magnitude ) 
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and phase is converted to degrees using the relation 
Phase = (180/7f) Phase 

Magnitude calculations for the single pole or zero can 
be written as 

Magnitude = [Realpart^ + ( w-Imaginarypart ]1/2 

and the phase calculation is 

Phase = Tan~^ [ ( w-Imaginarypart ) /Realpart] . 

For the whole system, the magnitude and phase are 
calculated for each pole or zero. Then the final 

magnitude is 

Magnitudegystem = Magni tudezeros/^^snitudepoigg 

and the phase is 

Phase gyg-{;gjjj = Phase^eros ~ P^^®®poles 

The plots of the magnitude and phase are shown on the 
sajne graph. Coordinate values of 0dB magnitude and 
-180° phase coincide in the graph. We know from 
control theory that phase margin is read at the zero 
crossover of the magnitude curve and the gain margin is 
read at the -180° crossover of the phase curve. These 
two values can be read directly from the graph. The 
Bode routine calculates the numbers required for the 
plots. The procedure "plot-Bode" converts the number 
to a graphical display. Within the Bode-plot routine 
is a call to the procedure "graph-menu". Graph-menu is 
called by all procedures which produce a graph. It 

provides a small menu offering the user the choice to 
add a title to the graph, print the graph on the 
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printer, print the number or quit and return to the 
menu . 

If the user selects the title to graph and it is 
completed, the plot is displayed with the title on the 
screen. The title block can be moved by using the 

cursor arrow keys ajid relocated anywhere on the screen. 
When the title box is moved where the user wants, the 
<enter> key must be pressed. Then the screen is frozen 
in position and the graph menu is recalled on the 
screen . 

The print numbers selection saves the current 
graphic screen and permits the user to print the 
numbers used to draw the graph. The numbers may be 
printed on the printer (this will use a lot of paper) 
or to a disk file. If the disk option is selected, the 
user caxL scan that file with a word processor or by 
using the DOS "type" command and examine the points of 
interest. This option is illustrated on example 3.1. 

EXAMPLE 3 . 1 

The example caji be stated as follows: Given plant 
transfer function 



Gp ( s ) = 



100 ( 0.02 S + 1 ) 



(s + 1 ) ( .1 s + 1 ) ( .01 s + 1 )- 



(5-1 ) 



was rearrajiged as a state variable equations. 



0 


1 


0 


0 


x(t) + 


0 


0 


0 


1 


0 


0 


0 


0 


0 


1 


0 


-1 00000 


-1 1 2000 


-1221 0 


-21 1 




1 



y(t) = [^10000000 200000 0 0^ x(t) 
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a) Obtain the Bode diagram of the above system. 

b) Mark the following on the Bode diagram, recording 
the numerical values. 

1 ) Gain crossover frequency 

2 ) Phase margin 

3) Phase crossover frequency 

4) Gain margin 

5) Resonant frequency. 

Solution : 

Given data entered to the program by the selecting 
of the "Input Plant Matrices" option in the 

input /Change menu. Then the program outputs can be 
seen from the Figure 3.3 and Figure 3.4. 




Figure 3.3 Open loop Bode plot for example 3.1. 
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Figure 3.4 Closed loop Bode plot for example 3.1. 



C. TIME RESPONSE PLOT 

It is desirable to see the system’s response in the 
time-domain to a typical input. Time response can be 
calculated and plotted with the SVS program as can the 
system input sinusoid, ramp, step or impulse. All 
these Inputs have user selectable amplitudes. Figure 
3.5 shows the time response parameter screen. 

The time response algorithm first converts the A, B 
and C matrices to the open loop transfer function, then 
into a discrete-time, state-space equivalent. The 
theory of the time response plot is not included here. 
Users who want more information about the subject 
should consult reference 4, or any other relevant 
textbooks . 

EXAMPLE 3.2 

For the given system in example 3.1, obtain the 
time response plot and mark the following on the time 
response plot, recording the numerical values. 
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1 ) Settling time 

2) Maximum overshoot for a step input. 

Solution : 

The program data has already entered for the Bode 
plot. The selecting "Time response Plot" in the 
graphics menu gives figure 3.6 as a problem solution. 



Time Response Plotting Parameters ****** 


What is your 


input to the system? STEP 


(S) S 




RAMP 


(R) 




SINE WAVE 


(W) 




IMPULSE 


(I) 


What is your 


input amplitude? 1 




Input one of 


these choices. Open (0) or 


Closed (C) C 


Input your simulation time to the system (99max) 1 .4 



Figure 3.5 Time response parameters screen. 



D. NYQUIST PLOT 

This section presents the Nyquist plot option. 
This selection gives open loop and unity-feedback 
closed loop Nyquist plot. The program first calculates 
open loop transfer function of the plant. Then the 
plot is obtained by calculating the magnitudes and 
phases angle of the transfer function for a specified 
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Figure 3.6 Time response plot for example 3-2. 



number of times with a specified increment of w. 
Before proceeding to the Nyquist plot, the program 
prompts the menu of Figure 3.7. This parameters screen 
allows the user to enter additional data for the 
plotting. The graphic window size is given 100X100 
scale if the user selects the big picture option. The 
select own size option asks the user to enter starting 
frequency, number of decades and X,Y coordinates 
maximum and minimum values for the plot. After getting 
the plot, the procedures are the same as with the Bode 
plot option. These are make a title to the graph, 
printer output and listing numbers (which are used to 
generate a graph) either to the printer or a specified 
file with a given drive naune . 
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MMM Nyquist Plotting Parameters 

Open (0) or (C) Closed loop plot? 0 

Graph window (B) or (S) Select your own size? S 

Input your first frequency to be plotted? . 1 

(Example; .01, 1, 100, etc. ) 

Input number of decades do you want plotted? 4 

X-Maximum 1 00 
X-Minimum -100 
Y-Maximum 1 00 
Y-Minimum -100 

Any changes to these parameters? ( Y / N ) N 
Press <F1 > to change previous entry 
Figure 3.7 Parameters selection for the Nyquist plot. 

EXAMPLE 3.3 

Figure 3.8 shows the Nyquist plot for the example 
3.1 A, B and C matrices. 




Figure 3.8 Example of the Nyquist plot. 
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E. ROOT LOCUS PLOT 



This option plots the loci of* the closed loop poles 
of a system with transfer function G( s ) = K N(s)/D(s) 
for varying gain. N(s) and D(s) are polynomials of the 
plant and the program calculates from the given A, B 
and C matrices. 

The user is prompted to enter the starting and 
ending gain values, maximum and minimum X,Y coordinate 
values and to select either positive or negative 
feedback. These plotting parameters are shown on 
Figure 5.9. After this input routine, the program 
assumes unity feedback and calculates root locations 
for varying gain and plots them. 



****** Root Locus Plotting Parameters 

Input STARTING value for the varying gain 0 
Input ENDING value for the varying gain 1 0 



X -Minimum 


-100 


X -Maximum 


1 00 


Y-Minimum 


-100 


Y-Maximum 


1 00 



Positive or Negative feedback? ( P / N ) N 
Any changes to these parameters? N 



Figure 5*9 Parameters screen for the root locus. 



EXAMPLE 3.4 

Figure 3.10 shows the Root locus plot for the given 
data on the example 3.1. 
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Figure 3.10 Example for the root locus plot. 




F. CHARACTERISTIC EQUATION ROOTS 

This option of the graphics menu allows the user to 
decide whether the system is stable or unstable by 
looking at the root location of the characteristic 
equation. The program gives again the unity-feedback 
closed loop characteristic equation roots of the plant. 
The illustrative program output of the example 3-1 can 
be seen in Figure 3.11. 

G. LOAD PLANT MATRICES FROM DISK FILE 

This last option of the graphics menu allows the 
user to get saved data from disk file directly instead 
of going via to the main menu and input/change menu. 
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Plant Characteristic Equation Roots ***** 
ROOTS OF THE NUMERATOR 

s[1]= -50.00 0.000 

ROOTS OF THE DENOMINATOR 

s[1]= -152.032 +3 0.000 s[2]= -69.170 +j 0.000 

s[3]= -4.899 +3 -32.893 s[4]= -4.899 +3 32.899 



Figure 3.11 Example of characteristic equation 

roots . 
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IV. SVS MAIN MENU 



A. INTRODUCTION TO THE SVS MAIN MENU 

In this chapter six options are presented which may 
be used for the analysis and design of control systems. 
These options are supported by the other options which 
are explained in Chapter II and III. 

The observability option is used to determine the 
observability index of the system. Another very 
similar program, controllability is used to determine 
the controllability of the system. The matrix 
mathematics option brings to the screen another menu 
selections. This program calculates an A matrix 
determinant, inverse, characteristic polynomial and 
eigenvalues. The last three options may be used to 
design optimal linear control systems. The pole 
placement option is useful in the design of linear 
control state variable feedback control systems. In 
the pole placement case, the control is computed by 
multiplying by a gain [K] , the difference between the 
reference input and a weighted (linear) sum of the 
state variables. The Luenberger observer design is 
used to design a combined observer-controller to 
achieve a given desired closed loop transfer function 
when some of the states are not accessible. 

Design of optimal control will minimize a given 
cost function which produces a scalar control . The 
prograun starts to work at the terminal time and works 
backwards in time. 



B. CONTROLLABILITY 

This option is used to determine the 
controllability of the linear time -invariant system. 
Consider the following continuous-time system 
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where 



x(t) = A x(t) + B u(t) 



(4-1 ) 



X = state vector 
A = plant matrix 
B = input matrix 
u = control input 



The system described by the above equation is said 
to be state controllable at a given initial time if it 
is possible to construct an unconstrained control 
signal which will transfer an initial state to any 
final state in a finite time interval [Ref. 5]. 

This requires an algebraic condition such that the 
rank r(C) of the controllability condition matrix 



C = [B! AB! . . . . ! A^"”'b] (4-2) 

is n, the order of the system. 



EXAMPLE 4 . 1 

Consider the matrices A and B, 





0 


1 


-1 






1 


0 




1 


0 


1 


9 


B = 


0 


1 




0 


1 


-1 






0 


0 



Determine if [A, 



B] is a controllable pair. 



Solution : Since A is 3x3 and B is 3x2, matrix C has to 
be 3x6. The program checked the rank of the 
controllability condition matrix. The result of the 
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program can be seen in Figure 4.1, the system is 
controllable . 





CONTROLLABILITY RESULT 


The Plant matrix A is : 


0 . 00000E+00 


1 .00000E+00 —1 .0000E+00 


1 . 00000E+00 


0.00000E+00 1 .00000E+00 


0 . 00000E+00 


1 .00000E+00 -1 .0000E+00 


The input Matrix B is : 


1 . 00000E+00 


0 . 00000E+00 


0 . 00000E+00 


1 . 00000E+00 


0 . 00000E+00 


0 . 00000E+00 


The system is 


controllable . 



Figure 4.1 Controllability program output. 



C. OBSERVABILITY 

In this section we determine observability index of 
the linear systems. Consider the unforced system 
described by the following equations: 

The (unforced) time invariant system 

X = Ax (4-3) 

with the observation vector 

y = Cx (4-4) 

The observability index, which is defined as the rank 
r(0) of the observability condition matrix 



C 

CA 

CA^-I 



is n, the order of the system. 



(4-5) 
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The program output for this selection is 
illustrated in Figure 4.2 with the A matrix as given on 
example 4.1 and the following C matrix. 



C = [ 1 0 0 ] 



(4-6) 





OBSERVABILITY RESULT 


The Plant matrix A is : 






0 . 00000E+00 


1 . 00000E+00 


-1 


. 0000E+00 


1 . 00000E+00 


0 . 00000E+00 


1 . 


00000E+00 


0 . 00000E+00 


1 . 00000E+00 


-1 


. 0000E+00 


The Output Matrix C is : 






1 .00000E+00 


0 • 00000E+00 


0 . 


00000E+00 


The system is 


observable with 


index 3. 



Figure 4.2 Observability program output. 



D. POLE PLACEMENT 

Any single-input single-output 
invariant system is described by 
equations : 



linear, time 
the following 



x(t) = A x(t) + B u(t) (4-7) 

y(t) = C x(t) (4-8) 



The plant is characterized as a block diagram in 
Figure 4.3. 
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Figure 4.3 Block diagram for the state variable 
representation . 



The closed loop nature of the system of Figure 4.4 
is showed by the presence of the controller. It 
generates the control signal u from the knowledge of 
the state variables. So, we can see that, except for 
the reference input [r] , the state of the plant x is 
the only information needed by the controller. The 



P PLANT ] 




Figure 4.4 General closed loop system with state 
variable feedback. 
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control input is oomputed by the multiplying by a gain 
K the difference between the reference input and a 
weighted (linear) sum of the state variables. As a 
mathematical expression, 

u= K [ r - ( k-|X-| + IC 2 X 2 + .... + )] (4-9) 

where the kj^ ’ s are referred to as feedback 
coefficients. The gain K is referred to as the 
controller gain. The equation (4-9) may be simplified 
by making it in the matrix notation 

u= K [ r - k'^ X ] (4-10) 

After all this notation, the graphical representation 
of the system configuration is shown in Figure 4.5. 




CONTROLLER 






Figure 4.5 Linear state variable feedback system. 
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We could obtain the closed loop transfer function 
Y(s)/R(s) from the state variables. Our approach is to 
force the system into the configuration shown in Figure 

4.6. 




Figure 4.6 Closed loop block diagram representation. 



We can see that Hgq(s) is given by 



Hg.(s)= 



kT x(s) 



' y(s) 

Since y(s)= C x(s), this expression becomes 

Heq(s)= f 

^ C x(s) 



(4-1 1 ) 



(4-12) 



Here it must be pointed out that the program 
calculates the gain K, so zero steady state error 
results from a step input. If the user wants other 
conditions, he may rescale K and k^ appropriately by 
hand. For example, assume he wants to desire to have 
the controller gain, K=K-| , but the program output shows 
that K=K 0 with the feedback coefficients k-| , k 2 , k^ . 
The procedure is then to modify the program outputs by 
setting K=K-| and setting 
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[ k-| k2 k^ ] 



( 4 - 13 ) 



kT = 



K 

K 



0 



This procedure does not change Y(s)/R(s) and satisfies 
the condition K=K-| . Under these notations and block 
diagrams, this option gives an open loop transfer 
function Y(s)/U(s) for the plant, the feedback transfer 
function Heq(s), the controller gain [K] and the 
feedback coefficients [k^] to achieve the desired 
closed loop characteristic polynomial. 

The desired closed loop characteristic polynomial 
is the denominator of Y(s)/R(s) and must agree with the 
order of the plant. The user can enter this polynomial 
either in coefficient form or factored form. If the 
user wants to enter the coefficient form, the 
coefficient of the highest degree term must be unity. 

The methodology for computing the result is: the 
coefficients of the denominator polynomial of 
Y(S)/R(s), which is the polynomial desired by the user, 
may be adjusted at will by proper selection of k and K. 
The closed loop zeros are equal to the open loop zeros. 
In other words, linear state variable feedback has no 
effect on the zeros of Y(s)/R(s). The program 
calculates the numerator of Heq(s) to achieve the 
desired characteristic polynomial. Note that the 
complete Heq(s) is calculated by taking the numerator 
of Heq(s) and dividing it by the numerator of 
Y(s)/U(s). In the program procedures, the coefficients 
of the characteristic polynomial are computed by the 
use of the Principle-Minor method. For the matrix 
inverse calculation, the program uses the 
diagonalization procedures. 



EXAMPLE 4 . 1 

The plant matrices of a third order system are 
given in the following state variable representation: 
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-1 . 0 


1 .0 


0 . 0 




0 


0.0 


0.0 


1 . 0 


x( t ) + 


0 


0.0 


-3.0 


0 . 0 




1 



y(t) = [ 1.0 1.0 0.0 ] x(t) (4-15) 



Find the feedforward (controller) gain [K] and the 
feedback coefficients [k^] required to achieve a closed 
loop transfer function of 



Y(s) 2 (s + 2) 

R(s) s^ + 4s^ + 6s + 4 



(4-16) 



The program result for this example is presented in 
Figure 4.7. 



POLE PLACEMENT RESULT 
The Plant matrix A is : 

-1 . 0000E+00 1 .00000E+00 0.00000E+00 

0 . 00000E+00 0.00000E+00 1.00000E+00 

0.00000E+00 -3.0000E+00 0.00000E+00 

The Input matrix B is : 

0 . 00000E+00 
0 . 00000E+00 
1 . 00000E+00 

The Output Matrix C is : 

1 .00000E+00 1 .00000E+00 0.00000E+00 

Denominator of Y(s)/U(s) - Descending powers of S : 

1.0000 1.0000 3.0000 3.0000 

Figure 4.7 Program output for the pole placement. 
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The poles of the Y(s)/U(s) are: 



REAL PART 




IMAGINARY PART 


- 1 . 0000 


+j 


0 . 0000 


0 . 0000 


+d 


-1 .7321 


0 . 0000 


+d 


1 .7321 



Numerator of Y(s)/U(s) - Descending powers of S : 

1 . 0000 2 . 0000 

The zeros of the Y(s)/U(s) are: 

REAL PART IMAGINARY PART 

- 2.0000 +3 0.0000 

Desired closed-loop Characteristic polynomial - 
Descending powers of S : 

1.0000 4.0000 6.0000 4.0000 

The roots of desired closed-loop characteristic 
polynomial are: 





REAL PART 




IMAGINARY PART 




-2 . 0000 


+d 


0 . 0000 




-1 . 0000 


+j 


-1 . 0000 




-1 . 0000 


+J 


1 . 0000 


Numerator 


of the Heq(s) 


is - 


Descending powe 


1 . 5000 


1 . 5000 


0.5000 




The roots 


of the Heq(s) 


are ; 






REAL PART 




IMAGINARY PART 




-0 .5000 


+d 


-0.2887 




-0 . 5000 


+J 


0.2887 



The feedback coefficients [k*^] are : 



0.5000 0.0000 1 .5000 

The gain [K] is : 2.0000 

Figure 4.7 Program output for the pole placement 
( continued ) 

The results shown in Figure 4.7 specify that the 
gain of the controller [K] is 2.0 and the feedback 
coefficients are k-|=0.5, k 2 = 0.0 and k^=1.5. 
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E. LUENBERGER OBSERVER DESIGN 

Consider a linear time invariant plant of the form 

x(t) = A x(t) + B u(t) (4-17) 

y( t ) = C x( t ) (4-18) 

Let a feedback control law for equations (4-17) and 
(4-18) will be 

u(t) = K [ r(t) - kT x(t)] (4-19) 



PLANT 




Figure 4.8 Luenberger observer block diagram. 



Assume without loss of generality that the plant is 
controllable and observable. Since the state of 
equation (4-17) is not directly available to implement 
equation (4-19), an observer of the form 

z(t)= F z(t) + G-j y(t) + G 2 u(t) (4-20) 

and replacing the true states with estimates yields 

u(t)= K [ r(t) - k^ z(t)] (4-21) 

k"^z(t)= h^ z(t) + q"^ y(t) (4-22) 

will be designed. The resulting closed loop system is 
shown in Figure 4.8. 
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Where 



X 

u 



y 

r 

z 



Gl 



and G 

q 



A 

B 

F 

2 

K 

T 



hT 



state vector 
control input 
output 

system forcing input 
estimated state vector 
plant matrix 
input matrix 

observer eigenvalues matrix 
observer gain matrices 
controller gain 

output feedback coefficients matrix 
observer feedback coefficients matrix 



Luenberger Observer Design Parameters 

Input degree of observer (10 max) 2 

Input the desired feedback coefficients in 

Factored <F> Form C 
or Coef f Icient . < C > Form 

Input observer characteristic polynomial in 

Factored <F> Form C 
or Coefficient <C> Form 

Press <ESC> to change it!. 

Then type your input with <ENTER> key 



Figure 4.9 Luenberger observer design parameters. 



The user has to enter the controller gain and the 
feedback coefficients which can be found by the use of 
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"pole placement" option. The prograjn asks the observer 
eigenvalues which are represented by the F matrix in 
the prograjn. The observer degree depends on the 
observability index. For exajnple, if the observability 
index r is the minimum integer , then the observer gain 
matrix [G] has order r. Simply the order of the 
observer, when the program prompted can be entered. As 
equal to or greater than (r-1). These input parameters 
are shown in Figure 4.9. 

EXAMPLE 4.2 

The example presented here for the fourth degree 
plant is taken from Desjardin’s [Ref. 2]. The plant is 
represented by the following equations. 

0 ] x(t) (4-23) 



x(t) + 



Solution : 

Stepi 

The x-| and X 2 are the only measurable states and we 
want to achieve following closed loop transfer 
function . 



0 

0 

0 

1 



u( t ) 



(4-24) 



y(t) = [ 10 20 0 



x( t ) = 



0 

0 

0 

0 



1 
0 
0 
1 5 



0 

1 

0 

-23 



0 

0 

1 

-9 



Y(s) 

R(s) 



1 

s^ + 6s^ + 17s^ + 28s +20 



(4-25) 
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The controller gain [K] and the feedback 
cosfficients required are found, as can be seen from 
Figure 4.10 by the use of the "Pole placement" option 
in the same menu. Results shown in Figure 4.10 
indicates that the feedback coefficients [k*^] are -3, 
-6, 13, 20 and the controller gain [K] equals unity. 



POLE PLACEMENT RESULT 



The Plant matrix A is 



0 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 



1 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 
-1 . 5000E+01 



0 . 00000E+00 
1 . 00000E+00 
0 . 00000E+00 
-2 .3000E+01 



0 . 00000E+00 
0 . 00000E+00 
1 . 00000E+00 
-9.0000E+00 



The Input matrix B is : 

0 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 
1 . 00000E+00 



The Output Matrix C is : 
2.00000E+01 1.00000E+01 



0.00000E+00 0.00000E+00 



Denominator of Y(s)/U(s) - Descending powers of S : 

1.0000 9.0000 23.0000 15.0000 0.0000 

Numerator of Y(s)/U(s) - Descending powers of S : 

20 . 0000 

Desired closed-loop Characteristic polynomial - 
Descending powers of S : 

1.0000 6.0000 17.0000 28.0000 20.0000 

Numerator of the Heq(s) is - Descending powers of S : 
-3.0000 -6.0000 13.0000 20.0000 

The feedback coefficients [ k ] are : 

20.0000 13.0000 -6.0000 -3.0000 

The gain [ K ] is : 1 .0000 



Figure 4.10 Pole placement result for the Luenberger 
observer design. 
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Step2 

The observability index is determined using the 
"Observability" option in the same menu. That results 
is shown in Figure 4.11. 



OBSERVABILITY RESULT 



The Plant matrix A is 



0 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 



1 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 
-1 . 5000E+01 



0 . 00000E+00 
1 . 00000E+00 
0 . 00000E+00 
-2 .3000E+01 



0 . 00000E+00 
0 . 00000E+00 
1 . 00000E+00 
-9 . 0000E+00 



The Output Matrix C is : 

I.00000E+00 0.00000E+00 0.00000E+00 

0.00000E+00 1 . 00000E+00 0.00000E+00 



0 . 00000E+00 
0 . 00000E+00 



The system (A,C) is observable with index 3. 



Figure 4.11 Observability result for Luenberger . 



Step3 

As can be seen from Figure 4.12, the system is 
completely controllable. 



CONTROLLABILITY RESULT 




The Plant matrix 


A is : 






0 . 00000E+00 1 . 
0.00000E+00 0. 
0.00000E+00 0. 
0.00000E+00 -1 


00000E+00 
00000E+00 
00000E+00 
. 5000E+0 1 


0 . 00000E+00 
1 . 00000E+00 
0 . 00000E+00 
-2 .3000E+01 


0 . 00000E+00 
0 . 00000E+00 
1 . 00000E+00 
-9 . 0000E+00 


The input Matrix 


B is : 






0 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 
1 . 00000E+00 








The system (A,B) 


is controllable. 





Figure 4.12 Controllability output for Luenberger. 
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taken from 
an observer 
For this 



Step4 

An observability index r=5 (result 
observability output) allows us to design 
of order equal to or greater than (r-1)=2. 



LUENBERGER OBSERVER RESULT 



The plant matrix A is : 

0.00000E+00 1.00000E+00 
0.00000E+00 0.00000E+00 
0.00000E+00 0.00000E+00 
0.00000E+00 -1.5000E+01 

The input matrix B is : 

0 . 00000E+00 
0 . 00000E+00 
0 . 00000E+00 
1 . 00000E+00 



0 . 00000E+00 
1 . 00000E+00 
0 . 00000E+00 
-2.3000E+01 



The output Matrix C is: 

1 .00000E+00 0.00000E+00 0.00000E+00 

0.00000E+00 1 .00000E+00 0.00000E+00 

The desired feedback coefficients are : 



2.00000E+01 
1 . 30000E+01 
-6 . 0000E+00 
-3.0000E+00 



0 • 00000E+00 
0 . 00000E+00 
1 . 00000E+00 
-9 . 0000E+00 



0 . 00000E+00 
0 . 00000E+00 



The Observer characteristic polynomial coefficients 
in ascending powers of S 

1 .40000E+01 7.50000E+00 1.00000E+00 

The F Matrix is : 



-7.5000E+00 1.00000E+00 

-1 .4000E+01 0.00000E+00 

The G1 Matrix is : 

8.55000E+01 2.92500E+01 

0.00000E+00 0.00000E+00 

The G2 Matrix is : 

-3 . 0000E+00 
-1 . 5000E+00 

The output feedback coefficients are : 
2400000E+01 8.50000E+00 

The compensator feedback coefficients are : 
1 .00000E+00 0.00000E+00 



Figure 4.13 The Luenberger observer design output. 
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example, we want to design observer eigenvalues of 
-3 . 5 and -4.0. 



Steo3 

The program results for all this input data is 
shown in Figure 4.13. From results, the complete 
system can be described as: 



0.0 1.0 0.0 0.0 




X-] ( t ) 




0 


0.0 0.0 1.0 0.0 




X2< t ) 


1 


0 


0.0 0.0 0.0 1.0 




X3( t ) 


"T 


0 


0.0 - 15 . - 23 . - 9 . 




X4( t ) 




1 



z( t ) = 



- 7.5 

- 14 . 



-1 .0 

0.0 



z^C t ) 
Z4( t ) 



+ 



85.5 

0.0 



29.25 




x-i (t ) 


0.0 




X2(t ) 



+ 



- 3.0 
-1 .5 



u(t ) 



u( t ) 



1.0 r( t ) - 
- [1.0 .0] 




F. DESIGN OF OPTIMAL CONTROL 

The system considered is described by the following 
state variable equation: 

x(t) = A x(t) + B u(t) (4-26) 

where 

X = state vector 
u = control input 
A = plant matrix 
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B = input matrix 



The design of optimal control will minimize the 
following cost function: 



and XT means transpose of X. 

The physical interpretation of J(N) is this: we 
wish to keep the state near zero without excessive 
control energy expenditure. The input parameters of the 
design of optimal control are entered in the beginning 
of the program. This screen can be seen in Figure 
4.14. After entering these parameters and the Q 
matrix, then the program calculates the feedback gain 
matrix which, when multiplied by the state vector, 
yields a scalar control. In the program procedures, 
the following equations were derived using dynamic 
programming, starting at the terminal time and working 
backwards . 



GT(k)=-[GAMMAT » P(k-1) * FI]/[GAMMAT * P(k-l) » 

GAMMA + R] , GT(0)=0 (4-30) 



N-1 



J(N) = 1/2[ XT(N) Q X(N) 




Q X(k) + R U2(k)] 
(4-27) 



k=0 



where the following are defined 



Q = noise covariance matrix 

N = time intervals over which the SUM is made 
R = scalar random input 



P(k) = PSIT(k) * P(k-1) * PSI(k) + Q + GT(k) 



R * G(k), P(0)=0 

PSI(k)=FI + GAMMA * GT(k), PSI(0) = 0 



(4-28) 

(4-29) 
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»»»**♦ Design of optimal control procedure 



Input number of time intervals for SUM procedure? 40 
What is your sample interval? 0-1 

What is the value of scalar R? 2.0 

For the following options which cost function 

do you want? 0 

COST=terminal+fuel+traj ectory or 
C0ST=terminal+tra3 ectory <0> 

COST=terminal+fuel or COST=terminal <1> 
where 

terminal= 1/2 XT(N) Q X(N) 

N-1 

trajectory= 1/2 ^ X(k ) Q X(k) 

k=0 

N-1 

fuel = 1 /2 u2(k) 

k=0 

Press <ESC> to change it!, 

Then type your input with <ENTER> key 



Figure 4.14 Optimal control parameters screen. 



For simplicity in programming, the following 
definitions are defined: 



Terminal= 1/2 XT(N) Q 


X(N) 


(4-31 ) 


N-1 






Trajectory= 1/2 ^ X ( k ) 


* Q * X(k) 


(4-32) 


If 







N-1 

Fuel = 1/2 5 R U2(k) (4-33) 

k=0 
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EXAMPLE 4.3 



Given the system equation ajid parameters described 
below find the discrete steady state gains for a sample 
interval 0.1 , scalar R is 1.0 and number of time 
intervals is 40. 

x( t ) = 

The graphic result of the program is shown in 
Figure 4.15 and Figure 4.16, the numerical output is in 
the Figure 4.17. 



0 

5 



1 

0 



x(t) + 



0 

1 



u( t ) 



(4-34) 




Figure 4.15 Optimal control graphic output . 



G. MATRIX MATHEMATICS MENU 

This option is used to get various calculations 
with the plant matrix A of a given linear control state 
variable system. The selection from the SVS main menu 
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This second menu 



brings the matrix mathematics menu, 
consist of following options: 

(1) The determinant of A matrix 

(2) The inverse of A matrix 

(3) The characteristic polynomial of A matrix 

(4) The Eigenvalues of A matrix 



Before selecting this option, the user must select 
''Input/Change Plant Matrices" option to enter the 
plant matrices. 




Figure 4.16 Optimal control graphic output ^2. 
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OPTIMAL CONTROL RESULT 

The order of the system is : 2 
The number of time intervals is :1.0E-01 
The scalar R is : 1.0000 
The sample interval is : 40 

The A matrix is: 

0.00000E+00 1.00000E+00 

5.00000E+00 0.00000E+00 

The B matrix is: 

0 . 00000E+00 
1 . 00000E+00 

The Q matrix is: 

1 . 00000E+00 0 . 00000E+00 

0 . 00000E+00 2.00000E+00 

The FI matrix is: 

1.02510E+00 1 .00835E-01 

5.04177E-01 1.02510E+00 

The GAMMA matrix is: 

5 . 02087E-03 
1 . 00835E-01 



MINIMIZATION OVER ALL STAGES 



( stages ) 


GT :2 


GT :2 


1 


-1 . 047E-01 


-2 . 031E-01 


2 


-3.243E-01 


-4.21 IE-01 


3 


-6.722E-01 


-6.665E-01 


4 


-1 . 1 64E+00 


-9 . 495E-01 


5 


-1 .809E+00 


-1 .275E+00 


6 


-2.597E+00 


-1 . 640E+00 


7 


-3 . 494E+00 


-2.031E+00 


8 


-4.437E+00 


-2.424E+00 


9 


-5 . 355E+00 


-2 . 794E+00 


1 0 


-6 . 1 85E+00 


-3 . 1 20E+00 


1 1 


-6 . 888E+00 


-3.390E+00 


1 2 


-7.452E+00 


-3 . O03E+00 


1 3 


-7.886E+00 


-3 .765E+00 


14 


-8.21 0E+00 


-3.883E+00 




-8.447E+00 


-3 . 969E+00 


1 6 


-8 . 61 7E+00 


-4 . 030E+00 


1 7 


-8.738E+00 


-4 .073E+00 


18 


-8.824E+00 


-4 . 1 04E+00 


1 9 


-8.884E+00 


-4 . 1 25E+00 


20 


-8. 927E+00 


-4 . 1 40E+00 


21 


-8. 957E+00 


-4 . 1 50E+00 


22 


-8 . 978E+00 


-4 . 1 57E+00 


23 


-8.993E+00 


-4. 1 63E+00 


24 


-9 . 004E+00 


-4 . 1 66E+00 


25 


-9.011 E+00 


-4 . 1 69E+00 


26 


—9 . 01 6E+00 


-4 . 1 70E+00 


27 


-9 . 020E+00 


-4. 172E+00 


28 


-9 . 023E+00 


-4 . 1 73E+00 


29 


-9 . 025E+00 


-4. 173E+00 



Figure 4.17 Optimal control numerical output. 
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V. CONCLUSION AND RECOMMENDATIONS 



The software, SVS, is written as a teaching 
learning tool for student use. It can be a nice tool 
for analysis and design of linear control state 
variable systems. The program is fully Interactive and 
menu driven. The user does not get lost in the 
program. The <Q> key always returns to the main menu 
or CTRL + C key interrupts the program. All options 
were tested solving several example problems. 

Hopefully, all "bugs" have been eliminated. 

Furthermore, the program can still be easily 
improved and expanded. These are listed below. 

(1) Adding the ability to handle the KALMAN 
filtering . 

(2) Add the capability of the program to handle 
discrete time systems. 
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APPENDIX A 



UTILITI FILES 

A set of Turbo Pascal input and menuing utilities 
was used widely in the program. These are public domain 
utilities and was copied from the LCS-CAD source 
diskettes. These utilities, called TURBO-UT . PAS were 
written by; 

Donald R. Ramsey 
Larry Romero 
727 Bunker Hill j^70 
Houston, Texas 77029 

ajad distributed through the public domain. The 
describing documentation are presented on reference 4. 
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APPENDIX B 



PROGRAM LISTING 

Appendix B is a listing files of the Turbo Pascal 
source code. These files make up the major modules of 
the SVS program. In general, most of the driver 
prograjns and including files are listed in the 
following pages. The source modules from the Borland 
International Turbo Graphix Toolbox and the 
input/output utility routines are not listed. 
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FILE: SVS.PAS 



Page 1 



Program Listing 

Program State Variable System( input , output ) ; 

(* *************** **** *^ * it **** ^'^ ^ * ^ ^ ^ Ik it It * *** it *) 

(* The following include files contain procedures *) 
(* to handle graphics call from the main procedures*) 
(* of this program. These include files are a part *) 
(* of Borland International’s Turbo Graphix Toolbox*) 
(* Which is a commercially available product. *) 

f * ********************************* *^ *********** * * i 

($I Typedef.SYS} (type & variable decleration) 

{$! Box. INC) {draw the box for the main menu) 

! * ********************************************** * 

* These are utility procedures and functions to * 

* help input. They are public domain programs. * 

* ********************************************** * 

($I Ut-Mod01 . INC ) (I/O procedures) 

{$I Ut-Mod02. INC) 

Procedure MainMenu; 
var 

I T 3 . Id i 

Help) File .InputFile ,MatrixFile , 

E lotf ile , ControlFile , ObserFile , 
uenbergFile , PoleFile , OptimalFile : 

Description : 

Okchoices : 



Integer ; 



File ; 
Str80: 
Set Ol* 



Char 



procedure ProgramExit; . 

(Displays warning about program end) 



be 



;in 
;lr scr 



Have you wanted ending the SVS program? ( Y/N ) 



Highvideo ; 
Msg (’Have yo 
S 1 0, 1 0) ; 
LowVideo ; 

Repeat 

option; ^ 

if not (Ch^in [’Y 
Until Ch in [ ’Y’ , ’N 
End; 






N’]) then Beep( 900 , 350 ) 



Procedure Menultem(pick : char : description: str80; 

color : integer ) ; , . 

{ allows easy selection of main menu colors) 

Begin 

Text Color ( color ^ \ 

WriteC’ ’:Tab,’(w; TextColor (white )j Write(pick); 
TextColor( color ) ; Writeln( ’ ) ’ , description ) ; 

End; 

Begin { Main Menu ) 

ClrScr: 

TextColor ( lightblue ) ; „ , x -t 

MsgC’ Naval Postgraduate School Ismail 

UNLU ’,7.24i; 

GoToXY( 21 , 4 ) ; TextColor (white ) ; 



Writeln( ’ 

Writeln( ’ 
Tab : =23 ; , 
Menultem( 
Wr i teln ; 
Menultem 
Menul tern 
Menultem 
Menultem 



) 



SYS MAIN MENU'*^»» ’)j 

{ show main menu ) 



Input / Change 




Plant Matrix Menu’, red); 

); 



’Graphics Menu ’ , 1 ightgrav ) 

’ Controllability’ ,yellow) ; 

’ Observability ’ .yellow) ; , 

’Luenberger Observer Design ’, yel low ) ; 
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FILE: SVS.PAS 



Prograim Listing 



Page 
.yellow) ; 



Menultem( ’ D Design of Optimal Control’ 

Menultem( ’ P ’ , ’Pole Placement ’ .yellow) ; 

Writeln ; 

Menultem( ’M’ . ’Matrix Mathematics Menu’ .lightgray); 
Writeln ; 

MenuItemC ’H’ . ’Help’ . lightmagenta) | 

Menultem( ’ Q ’ . ’ Quit the program’ . lightmagenta) ; 

Box; , ^ ^ , 

if not (blocki in [1..20] ) then 

OKchoices:= t’I’.;H%’Q’J 
else OKchoices :=[’!’ , ’Cl . ’0’ . ’D’ , ’M’ . ’P’ . 

’G’ . ’L’ . ^H’ . ^Q’] ; 

{sets legal choices depending on user) 



Repeat 

Option; ^ 

If not (ch in okchoices) then 
Begin 

Text Co lor ( red ) ; 

If blocki <= 0 then 
Begin ^ 

Beep( 900 . 550 ) ; 

Msg( ^WARNING: First 



End; 

End; Text Co lor ( white ) ; 
Until Ch in Okchoices; 
Case Ch of 

Begin 



INPUT .HELP or QUIT!’. 
1 ,25); 



’M’ 



5gin 

Asslgn( InputFile . ’ Input . com ’ ) ; 
Execute( InputFile ) ; 

End : 

Begin 

Assign(Matr ixFile . ’Matrix. com’ ) ; 
Execute(MatrixFile ) ; 

End : 

Begin 

Assign(Plotfile, ’plot.com’ ); 
Execute(plotf ile ) ; 

End j 
Begin 

Assign( ControlFi le . ’ Control . chn ’ ); 
ChainCControlFile ) ; 

End: 

Begin 

AssignC ObserFile , ’ Obser . chn ’ ) ; 
ChainC ObserFi le ) ; 

End: 

Begin 



AssignfLuenbergFile . ’ Luenberg . chn ’ ) ; 
Chain(LuenbergPile ) ; 



End, 
Begin 



H 



AssignC optimalf ile . ’ optimal . chn ’ ) ; 
Chain( opt imalf i le ) ; 

End: 

Begin 

Assign(Help1 File . ’ Helpl . chn ’ ) ; 
ChainCHelpIFile ) ; 

End ; 

Begin 

Assign( PoleFile . ’ Pole . com’ ) ; 
Executed PoleFile ) ; 

End ; 

Begin 

ProgramExit; 

if (ch=’Y’) or (ch = ’y’) then 

Exit := True; 

End; 
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FILE: SVS.PAS 



Program Listing 



Page 3 



End; {case} 

End; {main menu) 

* ***************************************** * 

* Program Execution part * 

egin {svs> 

ClrScr; {clear screen) 

Exit ;= False ;{ initialize boolean) 

Repeat 

if blockl <> 1 then 
Begin 

Blockl := 0; 

End; 

MainMenu; { call main menu until the user 
want to exit) 

Until Exit= True; 

Set Cap Num( ’ ’ , ’ ’ , ’ ’ ) ;Say Cap Num; 

“ , {set caps , insert , and num lock off) 

Blockl ;= 0; {reset problem) 

End. {svs) 
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FILE: INPUT. PAS 



Page 1 



Program Listing 

Program input ( input , output ) ; 

{* This program allows the user to enter , change , save 
and retrieve the problem for the whole options on the 
other menues. *) 



($I Typedef.SYS) (comman type & variable definitions) 

($I Ut-mod0 1 . INC ) {I/O utility routines) 

(SI Ut-mod02. INC) 

<$I Ut-mod03. INC) 



($I Box. INC) {draws menu box) 



var 

inputdatf ile , savef ile , retrievef ile 
changef ile ,help2file 



Procedure InputMenu; 
var i,Tab : 

Okchoices : 

Finished : 



Integer : 

Set Of Char 
Boolean ; 



f ile ; 



Procedure 



Begin 
TextCol 
Write ( 
TextCol 
End; 



MenuItemC Pick : Char ; Description : Str80 : 

Color : Integer ) ; 

{displays menu items in color ) 

or(color): ^ . s . ^ 

’ ;:tab, (’); TextColor( white ) : Write(pick) 
or (color); Writeln( ’ ) ’ , Description); 



Begin 

ClrScr: TextColor ( white ) ; Finished := False; 

GoToXY( 20 , 4) ; retriev: =true ; 

Write(^*** INPUT / CHANGE MENU *** ’); 

for i := 1 to 4 do writeln( ’ ’ ); 

Tab := 18: 

MenuItemC * I Input Plant Matrices ’,red); 

Menultem{ ’ C Change Current Plant Matrices yellow) ; 



Wr i teln ; 

Menultem( ’ S ’ , ’ Save Plant Matrices to Disk 
^ File ’ .yellow) ; 

MenuItemC ’ L Load Plant Matrices From Disk 

File ’ .yellow) ; 

Wr iteln ; 

MenuItemC ’ H ’ , ’ Help ’ , lightmagenta) ; 

MenuItemC ’ Q Quit to SVS Main Menu ’.lightmagenta); 
Box; TextColor ( white ) ; 



Set Cap num( ’ C ’ 
if block) <> 1 
else 



then OKchoices := 
OKchoices 





); 



Repeat {wait for user to input a keyprees) 
option ; ^ 

if not (oh in OKchoices) then 
begin 

beepC 900,350); TextColor ( red ) ; 
if blocki <> 1 then 

^®msg( ’WARNING: First INPUT , RETRIEVE . HELP or 

QUIT! ^,1.25); 

end; 

end; 

TextColor ( white ) ; 
until Ch in OKchoices; 
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FILE: INPUT. PAS 



Program Listing 



Page 2 



case ch of 
’ I ’ : begin 



’L’ 

’S’ 



Assign( inputdatf ile , ’ inputdat . chn ’ ) ; 

chain( inputdatf ile ) ; 
end ; 
begin 

Assign ( retr ievef ile , ’ retrieve . chn ’ ) ; 

Chain ( retr ievef ile ) ; 
end j 
begin 

AssignC savef ile , ’ save . chn ’ ) ; 

Chaln(saveflle) ; 
end; 



’ C ’ : Begin 
’H 



ssign( changef ile , ’ change . chn ’ ) ; 
Chain( changef ile ) ; 



End: 
Begin 



Assign(help2f ile , ’ help2 . chn ’ ) ; 
Chaln(help2flle^* 

End; 

’ Q ’ : begin 

Assign( svsFi le , ’ svs . com ’ ) ; 

Execute( svsFile ) ; 
r end: 

End; {case} 
end; 

{ Execution the input/change menu program. > 



Begin {main program) 
drive : = ’ C ’ ; 
repeat 

InputMenu ; 

Until Finished = True; 
End. {main program> 
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FILE: INPUTDAT.PAS Prograjn Listing 

Program input data( input , output ) ; 

{To allow the user to enter the A,B,C matrix from the 
this program) 

label 13; label 10 ; label 11 ; {label decleration for the 
GOTO statement) 

{$! Typedef . sys ) 



I 



:i Ut-mod01.inc 
;i Ut-mod02.inc 






{comman variable declerations ) 
{I/O routines) 



var Ans , Cont : Char ; 

Out : Text ; 

Stepping , Step , Steps , 

Tempi , Temp2 , i , 3 : Integer; 

Result : Real; 

inputfile : file; 

Begin {input data) 

1 3 : ClrScr ;~Wr i teln ; 



block 1 := 1 : textcolor (yellow) ; 

’Enter the degree of the plant: 



); 



size > 10 ) then 



{ initialize ) 



Wr ite( , 

Readlni size ) ; ^ 

if (size <= 0 ) or ( 
begin , 

beep{ 900 , 350 ) ; goto 13; 
end; 

for steps ;=1 to 10 do 

for stepping := 1 to 1 0 do 
A1 A[steps , stepping] ;= 0.0; 

Writeln( ’ Enter the elements of the A Matrix ’); 
wr iteln ; 

for steps : = 1 to size do 
Begin 

or stepping : = 1 to size do 



Begin 

Wr ite( ’ A 
Readln( A 
End;Writeln; 
End; 

Repeat 



( ’ ^steps, ’ , ’ .stepping, 
1 A [steps , stepping] ) ; 



) = ’ ) ; 



) ; wr iteln ; 
do 



v^lrScr; Writeln; 

Writeln(’The A Matrix is : 
for steps :=1 to size do 
Begin 

for stepping := 1 to size 
Begin 

Write( ’ A1A [steps , stepping] : 1 1 ) ; 

End; Writeln ; 

End: writeln; {prompt to any changes) 

Write( ’Do you want to change any element of’); 
writeC ’ the Matrix ? { Y / N T M; 

Read(Kbd, Ans ) ; writeln; {allows user to change 

, ^ , entered data) 

if { Ans =’Y’) or ( Ans =’y’) then 
Begin , X , X 

e:’); readln( i ) ; 
ange : ’ ) ; 

r< ■ ■ - - . - ^ . 

wr : 

A1 



writei ’Input the row to chang( 
write(’ Input the column to chj 
eadln( 1 ) ;wr iteln ; ^ ^ 

ite(’A(,i,’,’,3»’)= ’); readln( result ) ; 
A[i , 3 ] * =result ; 

’n’ ] ; 



End; _ 

Until Ans in [’N 



1 0 I CXrScr *wri“t©lri * 

{user inputs B mairix elements from the keyboard) 
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FILE: INPUTDAT . PAS Program Listing 

write( ’How many inputs do you have ? ’); 

readln(ni ) jwriteln; 
if ( ni < 1 ) or (ni > size) then 

begin 

beep( 900 , 350 ) ;goto 10; 
end; 

for steps := 1 to 1 0 do 
for stepping: =1 to 10 do 
. Brsteps , stepping] := 0 . 0 ; 

Writeln( ’Enter the elements of the B Matrix 



= 1 to size do 



Writeln; 
for steps 
Begin 

for stepping: =1 to ni do 
begin 

write( ; b(’ . steps, ’,’ ,stepping, ’ ) = ’); 

Readln( B L steps , stepping J); 
end;writeln; 

End; Writeln; 

Repeat 

ClrScr; Writeln: 

Writeln(’The B Matrix is: ’); writeln; 

for steps :=1 to size do 

Begin 

for stepping: =1 to ni do 
begin 

Write( ’ ’ ,B[steps , stepping] : 1 1 ) ; 

end; writeln; 

End; Writeln; 

WriteC ’Do you want to change any element of’); 
write?’ the B Matrix ? ( Y / N ) ’); 

Read(Kbd, Ans) ; 

writeln; {allows user to change B matrix element) 

if ( Ans =’Y’) or ( Ans =’y’) then 

Begin 

write? ’Input the row to change : ’); readln(i); 

write? ’Input the column to change ; ’); 

readln(i): writeln; 

writef’BC’, ’); readln? result ) ; 

B[i ] • =result ; 

End; . _ 

Until Ans in [’N’,’n’]; 

11:ClrScr; Writeln: ^ 

(user inputs output data from keyboard) 
write? ’How many outputs do you have ? ’); 

readln(no): . x 

if ( no <1 ) or ( no > size ) then 

begin ^ 

beepC 900 , 350 ) ;goto 11; 
end; 

for steps := 1 to 1 0 do 

for .stepping := 1 to 1 0 do 
Cl. steps , stepping] ; =0 . 0 ; 
wr* i "t G n * 

Writeln? ’Enter the elements of the C Matrix ’); 
Writeln ; 

for stepping := 1 to no do 
Begin 

for steps :=1 to size do 
Id Gsi n 

WriteC ’C( ’ , stepping, ’ , ’ , steps, ’ ) 

Readln( C [stepping , steps] ) ; 
end; writeln; 

End; writeln; 

Repeat 

ClrScr ; Writeln ; , 

WritelnC’The C matrix is : ’); writeln; 






) 
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Program Listing 



for stepping := 1 to no do 
Begin 

for^steps:=1 to size do 

^§rite( ’ ’ ,C[stepping, steps] : 1 1 ) ; 

end;writeln; 

End; Writeln; 



want . to change 



Wr i'fce 
write 

Read( Kbd, Ans ) ; writeln; 
if ( Ans =’Y’) or ( Ans =’y 



’ Do you 
’ the C 



Matrix ? ( Y / 



any 
N ) 



element 



of’ ) 



Begin 
wr ite( 



) then 



Input 



wr ite( ’ Input 



r eadl 
wcite 






the 

the 



writeln; 



row to 
column 



chnage : 
to change 



) ; , ^eadln( i ) 



writel ’ , 1 , ' , 
C » 3 J : =result ; 
End; 



3 > ’ )= ’); readln( result ) ; 



Until Ans in [’N’,’n’] 



Assign^ inputf ile , ’ input . com’ ) ; 

{ re-execute the input/change menu program} 
ExecuteC inputf ile ) ; 

End . { input data > 
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FILE: CHANGE. PAS Program Listing 
Program change data( input , output ) ; 



Page 1 



* This program allows the user to change A,B, C 

* matrices and their matrix order. 

abel lljlabel 12; label 14; 



<$I Typedef.sys) 



var 



Ut-mod01 . inc 
Ut-mod02 . inc 






{comman type & variable decleration) 
<1/0 utilities > 



Ans , temp , inputtype 
Stepping , Step , Steps , i , 
Result 



inputf ile 
result size 



Char ; 
Integer ; 
Real ; 
file; 
integer ; 



); 



Begin 

clrscr ;writeln; gotoxy( 1 , 22 ) ; 
invvideoC ’ Press <ESC> to chainge it!, 
gotoxyC 1 , 23 ) ; 

invvideoC ’ Then input your choice with <ENTER> key’); 
gotoxyf 1,2); 
textcolor ( lightblue ) ; 

writelnC ’ *** Change Current Plant Matrices 

Procedure 

TextColor(yellow) ; 

writelnC ’ = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = ’ ); 



to change 



MsgC ’Which matrix do you want 
textcolor C lightmagenta) ; 
msg( ’ 
msg( 
msgC ’ 

p 0'Q 0 3,'t 

input C ’A’ ,’A’ ,50,6j2,-^rue,F1 ,F10); 
temp : =copy< answec 
if not C^emp 



PLANT 

INPUT 

OUTPUT 




in [ ^ A ^ ’ 



until temp in [’A’ 
inputtype := temp; 

case^inputtype of 



B 



6’,’C’]) then beep<900,350); 
C’] ; 



A matrix is 
to size do 



’ ) ;writeln; 
do 



begin 

1 4 rClrScr ; 
writelnC ’ The 
for steps : = 1 
begin 

for stepping: =1 to size 
begin , ^ ^ , 

writeC ’ ’ ,A1A[steps, stepping] : 1 1 ) ; 

end;wr iteln; 
end;writeln; 

WriteC ’The order of the system is:’, 
size , ’ , ’ ) ; 

writeC’ Change ? <Y/N)’); 
repeat 

read(kbd,ch) _ 

until ch in I < Y ’ , ’ N ’ , ’ v ’ , ’ n ’ ] ; 
if (ch=’Y’) or (ch=’y’) then 
begin 

writeln; „ . . n 

writeC ’The order of the system is :’); 
readlnC result size); ^ 

if ( result_si^e< 1 ) or C result_size > 1 0 ) 
then ^11^^900^350). 14. 
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’ ) ;wr iteln ; 



end; 

size:=result size; 
end; 
repeat 

clrscr;wr iteln; 

Writeln(’The A Matrix is 
for steps :=1 to size do 
Begin 

for stepping := 1 to size do 

§rite(’ ’,A1A[steps, stepping] : 1 1 ) ; 
End; Writeln; 

EndtWriteln; 

Write( ’Do^y^u^w^t to change any element? 

Read(Kbg, Ans ); writeln ; 
if not CAns in [ ’N’,’n’]) then 

if ( Ans =’Y’) or ( Ans =’y’) then 
Begin 

write( ’Input row to change ; ’); 

readlni i ) ; 

writeC’ Input column to change : ’); 

readln(.i); writeln;, 
writef ’A( S i, d ; 

readlni result ) ; 

A1 A[i , d ] : = result ; 

End; ^ 

Until Ans in [’N’,’n’]; 
end; 

’ B ’ rbegin 

1 1 tClrScr ; 

writeln(’The B matrix is ;’);writeln; 

for steps:=1 to size do 

begin 

for stepping: =1 to ni do 
begin , _ _ , 

write( ’ B [steps , stepping] : 1 1 ) ; 
end; writeln ; 

©rid * wr* i "t © 1 n * 

Wri-te(’The number of input is : ’ ,ni 

Change ? (Y/N)’); 
repeat " ’ ’ ’ ’ 

until 
if (ch= 
begin 

writeln; {allow the user to change number 
, of input) 

writef ’The number of input is:’); 

readlni result size); 
if (result size<1) or (result size>10) 
then begin 

beep( 900 » 350 ) ; goto 11; 
end; 

ni:=result size; 
end; 
repeat 

clrscr; writeln; 

Writelnf’The B Matrix is : ’);writeln; 

for steps :=1 to size do 

Begin 

for stepping := 1 to ni do 
Begin 

Write( ’ ’, B [steps , stepping] : 1 1 ) ; 

End; Writeln; 

End;Wr iteln; 

(prompt the user for changes on B matrix) 



lange 7 (Y/N)’ ]; 

; read(kbd.ch) 
ch in [ ’ Y ’ , ’ N ’ , ’ Y ’ , ’ n ’ ] ; 
i=’Y’) or (ch=’y’) then 
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End; 

Until Ans in C’N’,’n’]; 



1 2 : ClrScr : 

writelnC’The C matrix is :’)jwriteln; 

for steps :=1 to no do 

begin 

for stepping: =1 to size do 
begin 

write( ’ C [steps , stepping] : 1 1 ) ; 
end; wr iteln ; 
end;wr iteln ; 

.Ilow the user to change number of outputs) 

’ »no, ’ ) ; 



wr*i"t©Xn * 

writet’^he Number of Output is :’); 
readlnC result size); 

if ( result_si'Ze< 1 ) or (result size>10) 
then begin 

beepC 900 , 350 ) ; goto 12; 
end; 

no:=result size; 



repeat 

clrscr ;wr iteln ;{ show C matrix elements) 
Writeln(’The C Matrix is : ’);writeln; 

for steps :=1 to no do 
Begin 

for stepping : = 1 to size do 
Begin , _ _ ^ 

Writei ’ C [steps , stepping] : 1 1 ) ; 

End; Writeln; 

End;Writeln; 

(prompt the user for ohanges on C matrix ) 




begin 



end; 




readln( result ) ; 
C [i » d ] : =result ; 
End; 



i » d ] : =result ; 
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end; 



Until Ans in [ ’ N 
end; 



Program Listing 
’n’] ; 



Assign( inputf ile , ’ input . com ’ ) ; 

ExecuteC inputf ile ) ; 

End . { change > 



{re-execute input 
program) 
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Program save_data; 

($I Typedef .SYS) 
($I Ut-mod01 . INC) 
{$! Ut-mod02 . inc ) 



var 

filename : str20; 
blockfile ttext; 
i , j : integer ; 

inputfile ‘.file; 
begin {save data) 
clrscr; HighVideo; 



{let user to change drive if necessary) 

Msg{’** Your Data Drive is . Press <esc> to change 
it ! , 1 0, 1 1 ) ; 

r* 0"p 0 

input ( ’ A ’ , copy( drive .1 ,1 ),32,11 ,2, true, FI ,F10); 
ch:= copy ( answer . 1 . 1 ) ; 

if not(cn in [ ’ A ^ ^ B ’ , ’ C ’ , ’ D ’ , ’ E ’ ] ) then 
beep{350 , 1 50 )i _ _ , 

until oh in [‘A’,’B’ ’ 



beep(^^y,.^_,T 

oh in C'A’ , ^6’ , ’C’ , »D’ , ’E’] ; 

Drive := concat ( oh , ’ : ’ ) ; 
clrscr; {prompt for filename to store data) ^ 

msg( ’ Input name of file to save data in File ’,1,10); 
wr* i G X n * wri'bGXn* 

wr iteln{ ’ * Your DATA disk must be in drive 

’ ,DrivG * 

input ( ’ A ’ , ’ ’ , 45 ’, 1 0 ,8 , true ,F1 ,F1 0 ) ; 

filename : = concat (Dr ive , copy( answer ,1 , 8 ) , ’ . svs ’ ) ; 



Assign( Blockfile , filename ) ; 
rewrite( Blockfile ) ; 

Wr itelnf Blockfile , size ) ; 
Writelni Blockfile ,ni ) ; 

Writeln( Blockfile ,no ) ; 
for i:= 1 to size do 

for .1: = 1 to size do ^ 

write (blockfile ,A1 A[i , j ] ) ; 
for i:=1 to size do 
for .1 ;= 1 to ni do 

write In (blockfile , B[i , j ] ) ; 
for i:=1 to no do 

for .1:=1 to size do^ 

write (blockfile ,C[i ,;) ] ) ; 



(Open file ) 



close ( Blockfile ) ; (close the file) 

Assign( inputf i le , ’input.com’ );{ re-execute 

input /change menu program) 
Execute ( inputfile ) ; 
end. (save data) 
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Program retrieve data; 

{$! Typedef.SYS) (variable decleratios) 

($I Ut-mod01 . INC > (I/O routines) 

($I Ut-mod02 . INC > 

($I Ut-mod03 . INC) 

($1 Dlrectry . INC ) (shows available data files) 



var 

readflle : text; (The text file user will use) 

filename : str20; 

readerror : boolean; 

llnecounter , i , i : integer: { A counter for the 

lines we read) 
inputf ile , PlotFl le : file; 

begin ( procedure retrieve data) 

Clr Scr ; HighVideo ; “ 

( Allow change of disk drive if desired) 
blockl j — 1 * 

Msg(’** Your data drive is . Press <ESC> to change 
it ! ** ’ , 1 0 , 1 1 ) ; 

?nput ( ’ A ’ , copy ( drive.1 ,1 ),32,11 ,2, true, FI ,F10); 



Ch:= copy ( answer , 1 , 1 
If notrcn in [’A^,’B 
beepC 900,3^021 






Until ch in 



B 



9 

’C’ 



C’,’D’,’E’]) then 



Drive := concat ( cA , ’ : ^ ; 



, ’D’ , ’E’ 



]; 



( Call directory to display eligible files) 
Dlrectory( drive , extension , filename , readerror ) ; 

If not ( reader ror ) then 
Begin 

Assignfreadf ile .filename ) ; 

( Open the file and read contents) 

Reset ( readf ile ) ; 

Llnecounter : =0 ; ^ (count the lines ) 

While not EOF( readf lie ) do 

(The built-in function EOF returns true) 
Beglndf the end of a file has been reached) 

Llnecounter : =linecounter+1 ;( Count the next line) 
Readln( readf lie , size ) : 

, ( Read itdnto variable line ) 

Readlnt readf ile ,ni 






Readln( readf ile ,no 
For i:=1 to size do 
For 1:= 1 to size do 

Read(readfile,A1A[i , J] ) ; 
For i:=1 to size do 
For 1:=1 to ni do 

Readln( readf ile , B [1 , j ] ) ; 
For i:= 1 to no do 
For j ; = 1 to size 
Read( readf ile . C | 

End; Close( readf ile ) ; 






d]); 



iy( 1500) ; 

. Wait fo 



End 
Else 
Begin 

Delj^., , 

/ait for directory error message and continue) 
Window ( 1 , 1 ,80 , 25 ) ; 

Clrscr ; 

End; 

if not (retriev) then 
begin 

AsslgnCplotf ile , ’Plot.COM’ ); 

Execute!, PlotFlle ) ; 
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end 

else 

begin ^ ^ , 

AssignC inputf ile , ’ input .com’ ) ; {re-execute , 

input/change menu program) 

ExecuteC inputf ile ) ; 

end: 

End. {retrieve data) 
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Program Plot: 

{ Program plot contains graphic programs. They are 
written by WOOD Roy, Jr. and modified for this 
program. > 

($I Typedef . SYS >{ type and variable decleration routine} 

($I Ut-mod0 1 . INC > {I/O routines } 

(|l Ut-mod02. INC) 

($I Ut-mod03 . INC) 

($I Box. INC) {drawing menu box ) 



var 

i 



■i 



help4f ile ,Retrievef ile ,Nyquistf ile , 



Integer ; 



timep lot file ,bodef ile ,rlocifile ,rootsf ile 



File; 



Procedure GraphicsMenu; 
var 

i , Tab 

Okchoice 

Finished 



Integer : 

Set Of Char 
Boolean ; 



Procedure Menultem( Pick:Char ;Description : Str80 ; 

Color : Integer ) ; 

{gives easy selection of input menu colors) 

Begin 

TextColor ( color ) : ^ , , 

Write (’ ;:tab,'(’)j TextColor (white ) : Write(pick); 

TextColor { color ); Writeln{’)’ .Description ) ; 

End; 



Be 



§lrScr ; TextColor ( white ); 
Finished := False : GoToXY( 20 , 4 ) ; 
Writelnf’ GRAPHICS MENU 



); 



{display graphics menu) 



Wr iteln ; 

Tab := 18 
Menultem( 

File^ 

Writeln ; 

Menultem( ’ C ’ , ’ Char acter istic Equation Roots 



’L’.’Load Plant Matrices From Disk 
»red) ; 



.yellow) ; 



’B’.’Bode Plot ’.yellow); 
’T’.’Time Response Plot ’ , > 

’N’.’Nyquist plot yel low ) ; ^ 

’R’.’Root Locus Plot ’ .yellow) ; 



Menultem 
Menultem 
Menultem 
Menultem 
wr* j.'boXn * 

Menultem( ’H’ , ’Help ’ . lightmagenta) ; 
Menultem( ’ Q ’ , ’ Qui t xo SVS Main Menu 



retr iev : =f alse ; 

{ read user choice from keyboard) 
L’ , ’B’ , ’T’ , ’N’ , ’H’ , ’Q’ , ’C’ ] ; 



, lightmagenta) ; 

.green;; 

Box; Writeln: 

TextColor (white ) ; 

GoToXY (40,22); 
repeat 
option; 

until oh in [ ’ R ’ , 
case oh of 
’R’ :begin 

assign( rlocifile, ’rloci. chn ’ ) ; 
chain ( rlocifile); 
end; 

’L ’ : Begin ^ . 

Assign( retr ievef i le , ’ retrieve . chn ’ ) ; 

Chain ( retr ievef i le ) ; 
end; 
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> rj\ i 

’B’ 
>N» 
’ C * 
’H’ 

End; 

End; 



begin 

Assign( timeplotf ile , ’ timeplot . chn ’ ) ; 
ChainC timeplotf ile ) ; 
end; 
begin 

Assign( Bodef ile , ’Bode. chn’ ); 
Chain(bodefile) ; 
end; 
begin 

Assignfnyquistf ile , ’ Nyquist . chn ’ ); 
ChainC Nyquistf lie ) ; 
end; 

Begin 

Assign( rootsf ile , ’ roots . chn ’ ) ; 

ChainC root sf ile ) ; 

End; 

Begin 

Assign(Help4File, ’Help4.CHN’ ); 
Chain(^Help4File) ; 

End; 

begin 

AssignC svsFile , ’ svs . com ’ ) ; 

ExecuteC svsFlle ) ; 
end; 



Begin (plot) 
drive : = ’ C ’ ; 

(initialize drive selection for load procedure) 

r* 0”p 0 

Gr aphicsMenu ; (repeadely call graphics menu until 

user selects to quit) 

Until Finished = True; 

End. (plot) 
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{ BODE PLOT is the driver program for the Bode 
plotting routines. It simply invokes Bode plot, when 
finished, returns back to the graphics menu. > 

Program Bode; 



$I 

SI 

SI 

SI 

SI 

$I 



Typedef . sys 
Graphix. sys 
Kernel . sys > 
Windows . sys 
Polygon.hgh 
Axis . hgh } 



1 

1 



{graphics routines) 



$I Ut-mod01.inc > {I/O routines) 
SI Ut-mod02.inc ) 

$I Ut-mod03.inc ) 



{$I Gr apMenu . inc ) (graph options menu) 

(SI PlotBode . inc ) {Bode plotting routine) 
{$I Boxuser.inc) 



type 



ary4 

aryo 



array [1..1i] of real; 
array [1..21J of integer; 



var 

plotfile tfile; 

($I Polynom. inc ) , 

{$I Rootf ind . inc ) 

{$I Bodeplot . inc ) 

begin 

BODEPLOT; {call 



{Polynomial routine) 

{Polynomial roots finder routine 

{Bode routine) 

the Bode calculation routine) 



) 



Assign(plotfile, ’PLOT.COM’ ); 

(re-execute the graphics menu routines) 
Execute (plotfile ) ; 



end. 
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Program time response: 

label 13; (label decleratlon for GOTO statement) 



Typedef .sys) 
Graphix . sys > 
Kernel .sys) 
Windows . sys ) 
Axis.hgh) 
Polygon .hgh) 



(graphics routines) 



$I Ut-mod01.inc 

§ 1 UT-mod02 . inc 
I Ut-mod03 . inc 



(I/O routines) 



($I GrapMenu. inc ) 



(graph option menu procedure) 



type 

ary4 

ary6 



array [1..11] of real; 
array [1**11J of integer; 



var 

Psi , Phi , A,Atemp 
temp, inputtype 
Offset , Slope , Tmax , 

RowSum , MaxRowSum , T , 

T1 , OldMaxRowSum , 

Plot time , Uinput , PhiX , 
hold, Ymax , Ymin , TPlot , 
Amplitude, Freq,y 
Factorial , Plot index , 

Nincr , code , i , 

J ' 3 » 3 n , kk , cnpoles , sizezero 

umpGraph , GoodNumber s , 
ClosedLoop , quit 
Cl , Xold, Xnext , Gamma 
GraphArr ay , Input array 
List 

numcoef f , dencoef f , cdencoef f 
kk 1 

plotf ile 



aryl s ; 
char ; 



real ; 



integer ; 

boolean ; 
ary3s ; 
plotarray ; 
text ; 
ary4 ; 
integer ; 
file; 



($I Polynom.lnc) (characteristic equation procedure) 
($I Boxuser.inc) 



Procedure Pr intGraphData; 

(this procedure dumps time-response data to printer) 
Begin 

LeaveGraphic ; Clrscr ; 
repeat 

TextcolorC white ): gotoxy(20, 1 0) ; 
writelnT’ PROGRAM OUTPUT OPTIONS 

gotoxy ( 20 ,13); 

writeln(’<P> Printer output^ 

Textcolor (yellow) ; gotoxy( 20,14); 
writeln(’ Check Your Printer! 

Textcolor (white ) ; gotoxy( 20 , 1 5 ) ; 
writeln(*<F> List to File name 
r(^20 , 1 6 ) ; 

,n( ’ on the current drive 



’ ); 
); 
); 

’ ); 
’ ); 



gotoxy 
writel-^, 

writeI^?*^<Q> 6uit ’) 

gotoxy^42 ,15) ; textcolor (yellow) ;wrlte( ’ "TIME.^ES" ’ ) 



gotoxv(28, 17) ; 

^If (ch=^F’;'or (ch=’f’) or (ch=’P’) or (ch=’p’) 
then 
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F> ) 



Program Listing 
or (oh = ’ f ’ ) then 



i'f (oh = 
begin 

gotoxy( 24,15); text col or ( red) ; 

writeC ’ PRINTING 

assign( list , ’ Time . RES ’ ) ; 
rewr ite( list ) ; 
end 
else 
begin 

gotox.yC 24,13) ; text col or ( red) ; 

writeC ^ PRINTING 

assign( list , ’LST : ’ ) ; 
rewr ite( list ) ; 
end 



); 



’ ); 



wr i"{;eln( list , ’ 

, Result 

wr itelni list ) ; 
writeln( list , ’ 



); 



TIME RESPONSE PLOT 



( input ) 



wr itelnf list , 
writeln( list ) ; 



, TIME (Sec) 
^ ( output ) 



R 



); 



do 

J ( 

f \ 



for i ;= 1 to 200 
wr iteln( list , ’ 

GraphArray 
InputArray [i , 2] : 1 

end; 

until oh in 
EnterGraphic ; 
swapscreen : 
close( list ) ; 
end ; 



IraphArray [i , 1 ] ; 1 0 
1.2J : 12 : 4 , 



5 ;’ 



Procedure Matrix Mult (Matrixl ,Matr ix 2 : aryl s ; 

var AnswerMatr ix: aryl s ; 
Order ; integer ) ; 



var i , J 



integer ; 



begin 

for i:=1 to order do 
for j;=1 to order do 

AnswerMatr ix[i , 3 ] ;= 0 ; 

(initialize the answer matrix) 

for i;= 1 to order do 
for 3:= 1 to order do 
for L : = 1 to prde;; do 
AnswerMatrix 



order do _ _ 

x[i..1J ;= AnswerMatr ix[i , 3 j 

Matrixl [i , L J *Matr ix2[L , 3 J ; 



end; 

Procedure Scalar Mult (Matrixl ; aryls; scalar 

var AnswerMatr ix : ary 1 s ; 
Order : integer ) ; 
var i,j : integer; 

begin 

for i:= 1 to order do 
for j:=1 to order do 

AnswerMatrix [ i , j J : = AnswerMatrix [ i , j J 

end; 



real ; 



* scalar; 



Procedure Matrix Vector Mult (Matrixl : aryls; 

Vector : ary 3 s; 

var AnswerVector : ary 3 s ; 

order : integer ) ; 
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var i ,3 : integer; 

begin 

for i:= 1 to order do 
begin 

nold:= 0; 
for 1 : = 1 to or 
hold:= hold + 

AnswerVector [ 

end; 
end; 



der do ^ ^ ^ , 

Matrixl [i , j ] ^Vector [j ] 
i J : = hold; 



Begin 

Initgraphic ; leavegraphic ; 

B o XU s 0 r * 

Characteristic Equation( A1 A , size jDencoef f ) ; 

{calculating denominator coefficient) 

for i:=1 to size do PSI [i , size] : =B [i , 1 ] ; 

for 33:=2 to size do 

begin 

for i:=1 to size do 
begin 

PSI ti , 3 J : =Dencoef f [kkl ] * B[i,l]; 

for 1;=’ to size^do ^ ^ ^ ^ , 

^ PSl[i,3]:= PSI[i,3] +A1A[i,l] * PSI[l.kk1]; 
end; 
end; 



for i;=1 to size do 

{calculating numerator coefficients > 

begin 

Numcoeff [i] : = 0.0; 

for j:=1 to size do ^ ^ ^ ^ ^ ^ 

Numcoef f [i J : = Numcoeff [i] + PSI[j,i] * C[l,i]; 

end; 

for i:=1 to size do {calculating numerator order > 
begin 

mi = siz04"1 * 

if Numcoeff [ni] <> 0.0 then goto 13; 
end; 

1 3 : sizezero ; =m-1 ; 

Clr scr ; 

Text Col or ( lightblue ) ; 

writelnC ’ ^ *** Time Response Plotting 

Parameters ***’); . . 

Text Co lor (yellow) ;writeln( ’==================’ 

MsgC’What is your input to the system? STEP (S) 

Msg( ’ ^6 amp (R) 

SINE wAv6 I'k) 

’ J 1 x8) 

Msg( ’ ■' " 



); 



Msg( ’ 



repeat , 



IMPl^L^E (l) 
’,1,9); 



Input ( ’A’ , ’S’ , 50,6 , 2 . true ,F1 ,F1 0 ) ; 
temp := copy( answer , 1 . 1 ) ; 
if not (temp in L ! ^ R ’ . ’ W ’ , ’ I ’ ] ) 
_Beep( 350 , 1 50 ) ; _ 

[’S’ ,^R’ , ^W’ , ‘i’] ; 



) then 



until temp in 



is input amplitude? ■ 

Input( ’N> , M ’ ,55,11 ,5, true, PI ,F10); 
val ( answer , Amplitude , code ) ; 



1 , 11 ); 
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case InputType of 
’ R ’ : begin 

Msg( ’What is your DC value? 



__ ’115) 

input (’’’N , ’0’"28' 15,5,true,F1 
val ( answer , Offset , code ) ; 
msg( ’What is your slope? ’, 1,15); 
Input ( ’N’ , ’1 ’ ,23, 15, 5, true, FI ,F10); 
val ( answer , slope , code ) ; 



); 



W' 



end j 
begin 
Msg( 



What is yotir frequency? 



end ; 



( rad/sec ) ’ , 1 , 1 3 ) ; 

Input ( ’N’ , ’ ’ ,35, 1^,5. true, FI ,F10); 
val ( answer , Freq , code ) ; 
end; 

Msg(’Open (0) or Closed (C) Loop Plot ? ’,1,17); 

P 0*0 0 3."t 

Input ( ’A’ , ’C’ ,40, 17,2 , true, FI ,F10); 
temp := copy ( answer , 1 , 1 ): 

if not (temp_in [’O’^’C’J) then beep( 350 , 1 50 ) ; 
until temp in [’0’,’C’]; 
if temp = ’C’ then ClosedLoop;= true 

else ClosedLoop;= false; 

Msg( ’ Input your simulation time for the system (99 

max r, 1,20); 

^®?nput( ’N’ , ’ ’ ,55, 20, 5, true, FI ,F10); 
val (answer ,Tmax , code ) ; 
if Tmax > 99 then beep( 350 , 1 50 ) ; 
until (Tmax <= 99) and ( Tmax >0) ; 

Boxuser ; 

if ClosedLoop then 
begin 

for i: = 1 to msixorder do CDenCoeffTiJ ;= 0.0; 

{ initialize) 

for i;=1 to^sizeZero + 1 do ^ 

CDenCoeffCiJ ;= Numcoeff[i] ; 

(C.L. denominator equals the sum of open loop 
denominator and O.L. numerator ) 



for i ; = 1 



do 



to SIZE +1 — 

CDenCoeffCiJ ;= CDenCoeff[i] + Dencoeff[i]; 
if size > sizezero then CNpoles:= size 

else CNPoles:= SizeZero; 

(Nploes shold always be greater, but to be safe) 



end 

else 

begin 

uNpoles : =Size ; 
for i;=1 to^size do ^ , 

CDenCoeffCiJ ; =Dencoeff [i] ; 

end; 

{ Calculation of new A matrix) 



for i : = 1 
begin 

for J :=1 
begin 

if J = 



end; 

end; 



to CNPoles-1 do 

to CN Poles do 

i+1 then A[i , i ] 
else A[i,j] 



1 . 0 

0 . 0 ; 
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{Calculation of new C matrix > 

for j:= 1 to CNPoles do 
begin 

A[CNPoles,o] := -CDenCoef f [^ ] ; 

3 n d 9 

for i:= 1 to CNPoles do 
begin 

if i > Siz 

if SizeZer 
Cl [i 
A 

end; 



eZero + 1 then C 



I 



1 [i] : = 0 . 0 

else Cl[iJ:= NumCoeff [i] ; 



= CNPoles then 



:= C1*fil + NumCoeff [SizeZero+1 1 * 
CNPoles,!]; 



{ Selection of sampling time interval ) 

Nincr : = 1000; 

T := ( TmcLx/Nincr ) ; 

A "t © mj^ J = A J 

Psi := A; {initialize psi to the value of the 
infinite series after the first) 



+ 1.0 

= T; 



Scalar_Mult( Psi ,T/ 2 , Psi , CNPoles ) 
for i:= 1 to CNPoles do {two terms I 
Psi [i ,i] :=Psi [i , ij ' “ 

Factorial := 2; T1 
repeat 
begin 

Factorial := Factorial * ( i + 1 ) ; T1 
Matr ix_Mul t ( A , ATemp , Phi , CNPoles ) ; 

{phi is used at temp) 

ATemp := Phi: {holding matrix to large array) 

Scalar_Mult ( Phi , ( T1 /Factorial ) , Phi , CNPoles ) ; 



+ A»T / 2 ! ) 



Oldmaxrowsum : = 0.0 



= T1 



for j:=1 to CNpoles do 
begin 

for m:= 1 to CNpoles do 
begin 

Psi [ j ,m] : =Psi [] ,m] + Phi[J,m]; 
end; 
end; 

Maxrowsum:= 0.0; 

{computes maLXrowsum as measure of change in last) 
for j:= 1 to CNpoles do{series term to be added.) 



begin 

rowsum:= 0.0; 

for m:=1 to CNPoles do _ _ 

rowsum:= rowsum + Psi [ j ,m] ; wr iteln ; 
if rowsum > maxrowsum then maxrowsum := rowsum; 



if {abs(maxrowsum-oldmaxrowsum)/maxrowsum)< 0.001 

then finished := false 
{quit when .1%change) else finished := true; 



o 1 dmaxr owsum 
end ; 

until Finished; 

Scalar Multf Psi , T , Psi 
Matr ix~Mult ( A , Psi , Phi 
for i^:= 1 to CNPoles. 

Phi[i,iJ := Phiri,i] 
for i:=_1 to CNPoIes 



= maxrowsum; 

CNPoles ) ; 
CNPoles ) ; 
do 

ao* 



~Gamma[i]”:= Psi t i , CNPolesI ; 
{single input system with B 



vector 



) 
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Plottime 
for i : = 

Ymax : 



< B= 



0 0 0 0 



0 1 



:= 0.0; Plotindex := 1; 

l[i J : = 



1 to CNf^oles do.XoldC 

{ init 

= 0.0; Ymin : = 0.0; 



prev . 



( transpose ) ) 
initialize ) 

. 0 ; 

state } 



Boxuser ; 



for N := 1 to Nincr do 

(begin calculating next state and y> 

begin 

(compute input at time Plottime) 
case Inputtype of 

’S’ : Uinput := Amplitude; 

’R’ : Uinput ;= Plottime * slope + offset; 

’I’ : if plottime = 0 then Uinput ;= amplitude 

else Uinput ;= 0.0; 

’W’ : Uinput := Amplitude * sin(freq * plottime) 

©nd * ^ c 3.S © )” 

Matf ix_Vector_Mult ( Phi , Xold, XNext , CNpoles ) ; 

(compute nev' states) 
for i : = 1 to CNpoles.do _ _ 

XnextLiJ:= Xnext[iJ + Gamma[i J ^Uinput ; 



y : = 0.0; 

for i : = 1 



begin 



to CNPoles do 



Tf’abs(y) < 1 . 0E07 then y:= y + Cl[i] 

' y:= 1.0E07; (max y limit 



else 

end; 

if SizeZero = CNPoles then 

y := y + NumCoef f [SizeZero+1 ] 
if y > Ymax then Ymax:= 



* Xnext^i] 



y ; 

if y < Ymin then Ymin:= y; 
if N mod 5=0 then (plot 
begin 

GraphArray .Plotindex, 1 . 
GraphArray .Plot index , 2. 
InputArray .Plotindex, 1 . 
InputArray [Plotindex, 2 
Plotindex := Plotindex + 
end ; 

Plottime := Plottime + T; 
end; 

Ymax := 1.1 * Ymax; 



* Uinput; 



every 5th point) 

= Plottime; 

= ?iottime; 

= Uinput; 

1 ; 

Xold ;= Xnext ; 



Initgraphic ; 

SelectWindow( 1 ) ; 

drawtext ( 20 , 20 , 1 , ’ 0 ’ ) ; 
drawtext (20,26,1 , ’ U ’ ) ; 
drawtext (20,32,1 , ’ T ’ ) ; 
drawtext(20,38, 1 , ’ P’ h 
drawtext (20,44,1 , ’ U ’ ) ; 
drawtext (20,50,1 , ’ T ’ ) ; 

drawtext ( 250 ,195,1,’ TIME( sec ) ’ ) ; 

NiceAxes ( 0 , tmaix .ymin .ymax , ’ ’ ) ; 

Selectwor ld(Wor IdNdxGlb ) ; SelectWindow( WindowNdxGlb ) 

DrawPolygon( GraphArray , 1 , -( Plot index-1 ) , 0 , 0 , 0 ) ; 
NiceAxes( 0,tmax, ymin, ymax, ’’ ); 

SetLineStyle( 1 ) ; (dashed line for input signal) 
DrawPolygon( InputArray , 1 , -( Plot index-1 ) , 0 , 0 , 0 ) ; 
SetLineStyleC 0 ) ; 
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repeat until keypressed; 
quit := false; 

P 0|3 0 

Graph Menu( ’ Time-Response ’ ,DumpGraph, quit ) ; , 

(calls prlnt/tixle menu) 
If DumpGraph then Pr IntGraphData; 
until quit; 

LeaveGraphic ; 

assign(plotfile , ’plot . com’ ) ; 

^ (re-execution graphics menu program } 
execute (plotfile ) ; 
end. 
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Program Nyquist; ^ 

label 1; (label decleration for goto statement) 



Typedef . sys > 
Graphix . sys > 
Kernel . sys 



(graphics routines) 



Windows . sys ) 
Polygon. hgh) 
Axis . hgh ) 

Ut-mod01 . inc 
Ut-mod02 . inc 
Ut-mod03 . inc 



Plotnyqs . inc 1 



$I Urap 
V>1 Plot 
$I Boxuser 



GrapMenu. inc 
. inc 
inc ) 



type 



ary4 

ary6 



($I Polynom. inc ) , 
($I Rootf ind . inc ) 



(I/O routines) 



(graph options menu) 
(Nyquist plotting routine) 



= array 
= ar r ay 



[1 . . 1 1 ] 

[ 1 . . 21 ] 



of 

of 



real ; 
integer ; 



(Polynomial routine) 

(Polynomial roots finder routine) 



var 



Code , I , Count .NumberDecades , 

StartDecade , EndDecade , one 

Wf , Wo , Wi , DeitaW. Gain 

PlotArrayl ,PlotArray2, 

MagPhaseArray.FreqArray 

ZMagn , ZPhase , PMagn , PPhase , Phase 

TempX , TempY 

temp 

OpenLoop , PicBig 

,kk1 ,m,l, sizezeros , CNpoles 
D( - ^ , 



integer ; 
Real; 

Plot Array ; 
real ; 
real ; 
char ; 
boolean ; 
integer ; 
ary4 ; 



)encoef f , Numcoef f , cdencoef f 
realpartpole , imagpartpole , realpar tzer b , 
crealpartpole , imagpartzero , cimagpartpole : ary3s ; 
PSI : aryls; 

plotfile :file; 

function Log( X : real ): real ; 

(computes the base-10 logarithm of X) 

Begin 

If X=0 then, 

Log := Ln(X 
End; 



)^og|-0.else 



10 ); 



function Expon( Y , X : real ): real ; 

(computes Y raised to X power) 

Begin 

Expon exp( X * (ln(Y))); 
end; 

begin 

Boxuser; one:=1; 

Characteristic equat ion( A1 A , size , Dencoef f ) ; 
for i:=l to si^e.do _ 

PSl[i,sizeI :=B[i,1] ; 
for jj:=2 to size do 
begin 

for i:=1 to size do 
begin 

1 : =size-j J + 1 ; 
kk 1 * = 1 “f* 1 * 

PSI ti , J J : =Dencoef f [kkl ] ^ B[i,l]; 
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* PSI[l,kkl]; 



* C[1 , 3 ] ; 



for 1 : = 1 
begin ^ 

|SI[i,j] :=PSl[i,j] +A1A[i,l] 
end; 
end; 
end; 

for i:=1 to size do 
begin ^ , 

Numcoeff [i] :=0.0; 
for j:=l to size do 
begin ^ , 

Numcoeff [i] : = Numcoeff [i] + PSI[J,i] 
end ; 
end; 

for i;=1 to size do 
begin 

m ; = siz0“f*1 i * 

if Numcoeff tm] <> 0.0 then goto 1 ; 
end ; 

1 : sizezer os : =m-1 ; 

Clr scr ; 

TextColor ( lightblue ) ; 

writeln( ’ *** Nyguist Plotting Parameters 

***' ) ; TextColor (yellow) ; 

writeln( ’ = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = ’ ); 
PicBig:= false; 



MsgC’Open (0) or (C) Closed Loop Plot ?’,5,5); 

’ ’,45,5,2 



repeat 

Input ( ’A 



rue ,F1 ,F1 0 ) ; 

-- then Beep(950,550) ; 

until temp in [’0’,’C’]: 

if temp = ’0’ then OpenLoop ;= true 

else OpenLoop ;= false; 

Msg( 



temp : =copy( answer , I , i ) ; 
if not(temp_in [’O’.’C’ 
il temp in [’0’,’C’J: 



Gr 



aph window (B) Big or (S) Select your own 

size? ’ ,5,7); 

r* 0TD 0 

, input ( ’ A ’ , ’ ’ , 60 , 7 , 2 , true ,F1 ,F1 0 ) ; , 

(sets flag OpenLpop if user ^selects the open loop) 

then beep( 550 , 1 50 ) ; 




else 



• cBig 
PicBig 



: = true 
:= false; 



if not (picbig) then 

^sg( ’ Input your^ f ir st ^ f re^i^ency to be 

Msg( ’( example . 01 , etc . ) ’ , 1 0 , 1 0 ) ; 

Input( ’N’ , ’ ’ ,50,9,8,true,F1 .F10) ; ^ ^ > 

Val (answer , Wo , code ); (Wo is the first plotted freq) 

Msg( ’ Input number of decades do you want 

plotted? ’ ,5,12); 

Input ( ’N’ , ’ ’ , 50 , 12, 2, true, FI ,F10); 

Val ( answer , NumberDecades , code ) ; 
end 
else 
begin 

wo : =0 . 001 ; 

NumberDecades : =8 ; 
end; 
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root_f Inder (sizezeros, Numcoeff , r ealp art zero , 
imagpar tzero , one ) : 

root_f inder ( size , Dencoef r , realpartpole , 
imagpartpole , one ) ; 

ain : =Numcoef f [sizezeros+1 ] ; 
or i:=1 to sizezeros+1 do 
begin 

Numcoef f [i ] : = Numcoeff [i] /gain; 

end ; 

clrscr ; Boxuser ; 

for i;=1 to maxorder do CDenCoeff[i] := 0.0; 

{ initialize > 

for i:=1 to^SizeZeros + 1 do, 

cpenCoeff[iJ := Numcoeff [i] * gain; 

(C.L. denominator equals the sum of open loop 
denominator and O.L numerator) 

for i:=1 to Size^+,1 do ^ ^ , 

CDenCoeff[i] := CDenCoeff[i] +Dencoeff [i] ; 



if Size > 



SizeZeros then Cnpoles : =Size 

else CNPOLES : =SizeZeros ; 



^ (compute new denominator roots) 

root_finder( Cnpoles ,CI)enCoeff .CRealPartPole , 
CImagPartPole , one ) ; 



StartDecade , 

.compute linear scale to plot ) 
EndDecade := StartDecade + NumberDecades ; 



trunc(Log(Wo ) ) ; 

(compute li: 

artDecade + Num 

(log numbers. Also figure step) 
Wf ; = Wo * Exponf 1 0 . 0 ,NumberDecades ) ; 

DeltaW := Expon( ( Wf /Wo ) , 0 . 01 25 ) ; 



Wi := Wo; 
for Count 



= 1 



to 81 do 

.(do 81 iter at ions ... arbitrary w) 

Begin 

if OpenLoop then 

(compute bode numbers if openloop and later if 
closed loop ) 

begin 

zMagn;=1.0; ZPhase:=0.0; 

PMagn : = 1 . 0 ; PPhase:=0.0; (initialize) 
for 1 : = 1 to SizeZeros do 
(compute magn and phase of zeros for freq step) 
tiegin , ^ r ^ 

ZMagn:=ZMagn * Sqrt(Sqr(RealPartZero[I] )+, 

Sqr(Wi-ImagPartZero[l] ) ) ; 



if RealPartZero[I] = 0.0 then 

ZPhase : =ZPhase+pi/2 . 0 else 
begin ^ , 

if realpartzero [i J > 0.0 then 

ZPhase := ZPhase - pi,+ ArcTan((Wi- ^ , 

ImagPartZero[i] )/( -Real Par tZero[i] ) ) 

else 

ZPhase : =ZPhase+ 



end ; 
end; 



ImagPartZero ealPartZero[l] ) ) ; 



for 1 := 1 to Size do 

(compute magn and phase of poles for freq step) 
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begin 



Magn:=PMagn * Sqrt ( Sqr ( RealPar t Pole [ I] ) + 

Sqr twi-ImagPart Pole L I J ) ) ; 



if* RealPartPole [i] = 0.(2) then 
PPhase:= PPhase+pi/2 . 0 else 
begin ^ , 

if* Real Part Pole [i] > 0.0 then 

PPhase : =PPhas§ -,bi + ArcTan(.(Wi- 



PPhase+^r^'^^^C 



RealPart Pole [I] ) ) 

Wi- ^ ^ 

ealPart Pole [i] ) ) ; 



ImagPartPole [i 

else 

PPhase:= 

ImagPartPole 

end; 
end; 

Phase := Frac((ZPhase - PPhase ) / ( 2*pi ) ) * (2*pi); 

, {Phase "modulo" 2Pi} 



,, , (Phase "modulo’’ 

TempX := abs( ( Gain»ZMagn/PMagn ) *cos( Phase ) 
TempY := abs( ( Gain*ZMagn/PMagn ) *sin(_ Phase ) 
^ and {TempX > 100j '' 



if (PicBi 
if (PicBig 



abj 

if) 



and (TempY > 100 



then TempX 
then TempY 



: = 1 0(3 ; 
: = 100 ; 



If Phase<0 then Phase := Phase+( 2*pi ) : 

If (( Phase > (pi/2) ) and ( Phase<13*pl/ 2 ) ) ) then 
^ TemnX;= -TempX; 

If ((Phase>pi) and ( Phase< ( 2*pi ) ) ) then 

_ _ TempY ;= -TempY; 

MagPhaseArr ay f count , 2 J "" * 



MagPhaseArray [count 
PlotArrayl tCount , 1 J 
PlotArray2 [count » 1 J 

PlotArrayl [Count, 2] 
PlotArray2 [count , 2] 

FreqArr ay [Count , 1 ] 



• — Ph3,S© J 
:= Gain*ZMagn/Pmagn ; 
TempX ; 

TempX ; 

TempY ; 

-TempY ; 



= wi 



Wi 

end 



= Wi * DeltaW; 



{increment freq step) 



©Is© 

{perform same 

begin 

ZMagn : = 1 
for i : = 
begin 



steps as above if closed loop requested) 

0 ; ZPhase ; =0 . 0 ; PMagn ; = 1 . 0 ; PPhase : =0 . 0 ; 

1 to SizeZeros do 



Magn : =ZMagn * Sqrt(Sgr (RealPartZero[Il^+^ 

■- -Sqr ( Wi-ImagPartZero [IJ ) ) ; 
if RealPartZero[l] =0.0 then 
ZPhase : =ZPhase+pi/2 . 0 else 
begin ^ ^ 

if RealPartZero[I] > 0.0 then 

ZPhase : =ZPhase - pi, + ArcTan( ( Wi- _ 

ImagPartZero[i] ) / ( -Real Part Zero [i ] ) ) 

else 



end; 

end; 

for i 
begin 



ZPhase : =ZPhase+ArcTan{ (Wi- 
ImagPartZero[i J) / ( -P 



Real Par tZero [i] ) ) ; 



= 1 to CNpoles do 



Magn:=PMagn * Sqrt(Sqr(CRealPartPoleril)+ 

Sgr rwi-CImagPartPole[IJ ) ) ; 
if CRealPartPole [I J = .0.0 



Wi-ClmagPartPole[I 

= 0.0 then 

PPhase : =PPhase+pi/2 . 0 else 
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begin 

if 

Clma^ 



CRealPartPole[I] > 0.0 then 
PPhase : =PPhase - pi + ArcT 
Part Pole [ij)/(-CRealPart Pol 
se 

PPhase : =PPhase+ArcTan( (Wi- ^ 

ClmagPartPoleliJT/ ( -CReal Part Pole [i] ) ) 






end; 

end; 



Phase := Frac((ZPhase - PPhase )/( 2*pi ) ) * (2*pi); 

(Phase ''modulo” 2Pi> 

TempX := abs( (Gain*ZMagn/PMagn)*cos( Phase) ) ; 

TempY := abs( ( Gain*ZMagn/PMagn ) *sin( Phase )) ; 
if (PicBig) and (TempX > 100) xhen 'TempX := 100; 
if (PicBig) and (TempY > 100) then TempY := 100; 



If Phase<0 then Phase := Phase+( 2*pi ); ^ ^ ^ 

If ( (^Phase> (pi/2) ) and ( Phase< ( 3*pi/^ ) ) ) then 

)X;= -Tempy ; 
len 

:= -TempY; 



If ((Phase>pi) and ( Phase<2*pi ) ) tS« 

•TempY 



PlotArrayl [Count , 1 
Plotarray2 [count , 1 

PlotArrayl [Count, 2 
PlotArr ay2 [count , 2 



] 

3 



TempX ; 

TempX ; 

TempY ; 

-1.0 * TempY ; 



Wi ;= Wi * DeltaW; 
end; 

Plot_Nyquist ( Start Decade , EndDecade , 

NumberDecades , Fr eqAr ray , PlotArrayl . 
Plotarr ay2 ,MagPhaseArr ay , PicBig , OpenLoop ) ; 



Assign( Plot file, 'Plot. COM' ) ; 
Execute(plotf ile ) ; 

end. (Nyquist) 
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Program Root_Locus; 

{ This program plots the root loous of the plant 
label !;<. label decler at ions for the goto statement) 



$I Typedef . sys > 
Graphix.sys) 
Kernel . sys ) ^ 
Windows . sys ) 
Polygon. hgh) 
Axis .hgh) 

{$! Ut-mod01.ino ) 
\$I Ut-mod02 . inc ) 
{$! Ut-mod03.inc ) 

<$I GrapMenu . inc ) 

type 



{graphics routines) 



{I/O routines) 



{graph options menu) 



Var 



n . .11] 
[ 1 . . 21 ] 



ary 4 = array 

ary 6 = array 

I , J , code ,LineCount 
Pi ox Pole, Plot Zero 



of real ; 
of integer ; 



integer ; 
PlotArray ; 



PlotRealPole , PlotImagPole 
DeltaGain, Star t Gain , EndGain , 

Variable Gain , Xmin , Xmax , Ymin , Ymax , gain : 
Neg_Feedback , dumpgraph , quit 

list 

1 j ,kkl ,m, 1 , size zeros , one 
holdpoly , dencoeff ,numcoeff 
realpartpole , imagpartpole , realpartzero , 
imagpartzero 
PS I 

plotf ile 

Polynom . inc ) 

Rootf ind. inc ) 



ary5s ; 

Real ; 
Boolean ; 

text ; 
integer ; 
ary4 ; 

ary3s ; 
aryl s ; 
file ; 



Procedure Pr intGraphData; , 

{dumps root locus data to printer) 

Begin 

LeaveGr aphic ; Clrscr; 
repeat 

Text CO lor {white ) ; 
gotoxy( 20,10); 
wr iteln{ ’ * * * 

gotoxy(20 » 13 ) 
writelnf 



PROGRAM OUTPUT OPTIONS * 



); 



’ ); 



_ _ *<P> ^’r inter output 

Text CO lor (yellow) ; gotoxy( 20,14); 
writeln( ’ Check Your Printer! 

Textcolor (white) ; gotoxy( 20 , 1 5 ) ; 
writein(’<F> List to File name 
gotoxy(20 , 1 6 ) ; * ^ . 

writelnC’ on the current drive 

gotoxy(20, 17) ; 

writeln{ ’ <Q> 6uit ^ ^ 

gotoxY{42, 15) ; textcolor(yellow) ; 
writeC ^ "Rioci . RES" * ) ; gotoxy( 28 , 1 7 ) ; 
readC kbd , oh ) : , n / , 

If (oh = *fM or (oh = ’f’) or (oh = ’P’).or 

( oh = ’p * 

b ©sin 

if (oh = ’F’) or (oh = ’f’) then 

begin 



’ ); 
’ ); 
’ ); 



) then 
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gotoxyf 24- , 1 5 ) ; text col or ( red) ; 

writeC ^ PRINTING 

assign( list, ’Rloci. RES ’ ) ; 
rewr ite( list ) ; 
end 
else 
begin 



f ^ PRiNT^^rG 



gotox 

write ^ , 

assign( list , ’ LST : ’ ) 5 
rewrite(llst ) ; 
end; 



; text col or (red) ; 



Page 2 
’ ) ; 



); 



LineCount := 0; 
wr itelni list ) : 
writelnt list . ’ ZEROS '); 
wr iteln( list ) ; 
wr ite( list , ’ 
wr iteln( list , ’ REAL 

^ IMAGINARY ’ 1 ; 
wr iteln( list ) ; LineCount 



’ ); 

= LineCount + 6; 



for i := 1 to sizeZeros 
begin 

writeln( list , ’ 



do 

’ , Real Par tZer o [i] : 1 0 : 3 , 

* , ImagPartZero [l] : 1 0 : 3 ) ; 



LineCount := 
end: 

wr itelni list): 
wr itelni list . ’ 
wr lteln( list ) ; 
write! list , ’ 
wr iteln( list , ’ 

wr iteln( list ) ; 
Variable Gain 



LineCount + 1 ; 

writeln( list ) ; 
POLES ’ ) ; 



GAIN ’ ) ; 

REAL 

IMAGINARY ’ ) ; 
LineCount := LineCount 
:= StartGaln; 



+ 7; 



(compute root locations for varying values of 
gain and print them) 

DeltaGain := ( EndGain-StartGaln ) /50 ; 

For J : = 1 to 50 do 
Begin 

HoldPoly := dencoeff; 



If Neg_Feedback then 
For I:= 1 to, sizeZeros+1 do, 
HoldPoly LI] := HoldPoly [I] + 
( gain*Var lable Gain * nnm 
else “ 

For I : = 



numcoeff [I] ) 



1 to sizeZeros +1 do^ , 

HoldPoly LIJ := HoldPoly[I] +( gain * , 

Var lable_Gain * numcoef f [ I] ) ; 

Root_Finder( size .HoldPoly , PlotRealPole , 

^ Plot ImagPole , one ) ; 
wrlteln( list, Var table Gain: 10:4); 



LineCount := LineCount + 1 
for 1 := 1 to size do 

begin 

wr lteln( list , ‘ 



, 1:2 



Plot Real Pole[i] Ll0:!^ 



, Plot ImagPole 



tl?: l6: 



LineCount := LineCount + 1 
end : 

wr itelnC list ) ;LineCount := LineCount + 
if LineCount > 50 then 



3); 
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begin 

wr iteln( list , chr(12)); 

LineCount : = 0 ; 
end j 

var iable_gain : =variable gain+deltaGain ; 
end ; ~ 

end; 

until ch in [’Q’,’q’]; 

EnterGraphic ; 
swapscreen : 
close( list ) ; 
end; 



Begin 
ope; 



. 1 . 
. 2 . 
.3. 
A. 



P .5. 



. 6 . 

7 




Root Locus Input handler driver} 
5506N01 001 -01 01 01 • 



5508N01 002-01 01 03 ’ 
1 51 2N00503-01 01 01 ’ 
151 3N00504-01 01 03 ’ 
1 51 4N00505-01 01 01 ’ 
1 51 5N00506-01 01 03 ’ 
451 7A00207T01 01 01 ’ 



Clrscr; TextColor ( lightblue ) ; 

writeln( ’ ROOT LOCUS PLOTTING PARAMETERS 

TextCol^r(Yellow) ; 

wr iteln ; wr iteln : wr iteln; 

TextColor (green ); 

writelnC ’ Input STARTING value for the variable 
gain: ’ ) ; 

wpi“tsXn * 

writeln{ ’ Input ENDING value for the variable 
eain: ’ ) ; 

TextColor (yellow) ; 
writeln; writeln; 



); 

); 



writeln: TextColor (green ) ; 

writelnC ’ X-Mlnimum: 

writelnC ’ X-Maximum: 

writelnC ’ Y-Minimum: 

writelnC ’Y-Maiximum: ’); writeln; 

writelnC ’ Positive or Negative Feedback? (P or N):’); 



Input Handler (’ N01 07 Escape ) ; 

~ (prompts for NEW inputs) 

wr iteln ; writeln ; 

■ 3 parameters? 

^(Y or N): ’ 

,ui.iuo,Ai,Fl0): 

then Input Handler (’ C01 07 ’ .Escape ) ; 
(prompts for changes) 



W X X u 6 X XX • Vr X V 6 X XX f 

writeln( ’Any changes to these 

Input ( ’A’ , ’ ’ ,45,19.2^true,F1 
If answer=’Y’ ' ” 



Val( f ilvar [ 1 ] , St art Gain , code ) 
(converts input strings 
EndGain.cQde) ; 
Xmin , code 
Xmax , code 
Ymin , code 
Ymax , code 



Val 

Val 

Val 

Val 

Val 



filvar .2 
fllvar .3 
filvar [4 
filvar 
filvar 






into numeric values) 



If copy ( filvar [7] , 1 

else 



1 ) < > ’ N ’ then 

Neg feedback := 
Neg feedback := 



false 
true ; 



INITGRAPHIC :( define values for graphics routine) 
N I CEAXESCxmin, xmax, ymin, ymax, ’ ’ ); 
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Character i St ic_equation( A1 A , size , dencoef f ) ; 
(calculate dencoeff from the plant matrices) 
for i:=1 to size do 

(calculate numcoeff from the given plant matrices) 
begin 

?Sl[i,size] :=B[i,l] ; 
end; 

for jj:=2 to size do 
begin 

for i:=1 to size do 
begin 

T:=size-dd+1 ; 

Iclc 1 * — H 1 * 

PSI ti » d J : =dencoeff [kkl ] * B[i,1]; 

for l:=i to size do 

°®fsi[i,d] :=PSI[i,d] +A1A[i,l] » PSI[l,kkl]; 
end; 
end; 
end; 

for i:=1 to size do 
begin 

numcoef f [i J : =0 . 0 ; 
for j:=1 to size do 
begin 

numcoeff [i] : =numcoeff[i] + PSI[j,i] * C[1,j]; 
end; 
end; 

for i:=1 to size do 
begin 

m * *“i * 

if numcoefftm] <> 0*0 then goto 1; 
end ; 

1 :sizezeros:=m-1 ; 

(calculate the zeros of the system) 
root_f inder ( sizezeros , numcoeff , real part zero , 
imagpartzero , one ) ; 

gain : =numcoef f [sizezeros+1 ] ; 

(convert highest degree numerator coefficient into 1 

for i:=1 to sizezeros+1 do_ 

numcoeff [ij : =numcoeff [ij/gain; 

For I:=1 to sizeZeros do 
Begin 

PlotZer 



•oCl, 1 ] 
■QU.2] 



PlotZer 
end; (for 
Case sizeZeros of 
0 



:= Real Par tZer o [ I] ; 
:= ImagPartZero[IJ ; 



1 



begin 

PlotZero .2 , 1 . : = PlotZero 
PlotZer o .2 , 2. :=PlotZero 
PlotZer o .3 , 1 . : = PlotZero 
PlotZero (3 , 2 J :=PlotZero 



1 



2 : 



1»2 ; 

NICEAXES(xmin,xmax jymin.ymax’ ’ ’ ); 
DrawPolygonC PlotZero , 1 , -3 , -3 , 3 , 0 ) ; 
end; 
begin 



Pl( 



•otZero [3 , 1 ] : = 
PlotZero [3 , 2 J : = 



=PlotZero [1 , 1 1 ; 
PlotZeroL 1 » 2 J ; 



NICEAXES(xmin.xmax,ymin,ymax, ’ ’ ) : 

Dr awPolygon( PlotZero ,1 , -3 , -3 , 3 , 0 ) ; 
d; 

NICEAXES( xmin , xraax ,ymin , ymax , ’ ’ ); 
DrawPolygon( PlotZero , 1 , sizeZeros ,-3,3,0) 



end; 
else 
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end; {case) 

Variable Gain := StartGain; 

DeltaGain := (EndGain-StartGain)/50 ; 

{divide gain to plot 50 points) 



For J : = 1 to 50 do 

{calculate and plot 
Begin 

HoldPoly := dencoeff; 

If Neg Feedback then 
For I:= 1 to^ SizeZeros+1 do, 
HoldPolyfl] ;= HoldPolyri] 



50 points per graph) 



(gain*Var iable_Gain * numcoeff [I] ) 



else 
For I : = 
HoldPol 




+ 

am * 



numcoeff [ I] ) ; 



Root_Finder ( size , HoldPoly , PlotRealPole , 
Plot ImagPole , one ) ; 



For 

Be 



I:=1 to size do 

{fill plotting matrix with poles) 
in 



[hi] 



’lotPole 
PlotPole 
end; 

AxisGlb := true; 

Case size of 

{artificially fill 
than 3 points) 



:= PlotRealPole 
:= PlotImagPole 



[hi 



plotting array if fewer 



0 

1 



be^in 



ot Pole 
PlotPole 
PlotPole 
PlotPole 



. 2 , 1 . 

. 2 , 2 . 

.3,1. 

3,2. 



= PlotPole 
=PlotPole 
=Plot Pole 
=PlotPole 



1 , 1 
1 ,2 
1 , 1 
1 ,2 



NICEAXES(xmin.xmax,ymin,ymax, ’ ’); 
DrawPolygonC PlotPole ,1 ,-3,-1 ,3,0); 
end; 
begin 

^lotPoleC3,1]:=PlotPole[1 ,1]; 

“,2J : =PlotPole[l ,2] ; 



PlotPole r 3 



end; 

else 



NICEAXES(xmin.xmax,ymin,ymax, ’ ’ ); 
DrawPolygon{ PlotPole , 1 ,-3,-1 ,3,0); 



NICEAXES(xmin,xmax,ymin,ymax, ’ ’ ); 
DrawPolygon( PlotPole , 1 , size , -1 ,3,0 



©ncL* vC^s©/ 
var iable_gain : =var iable_gain + Deltagain; 
end; 

Repeat until KeyPressed; 

quit := false; 
repeat 



Graph Menu ( ’Root Locus ’ .DumpGraph , quit ) ; 
“ alls print/title menu) 

then Pr intGraphData; 



If DumpGraph 
until quit; 

Leave Graph 1 C ; 
AssignCplotf lie , ’plot 
Execute{plotf ile ) ; 
end. {root locus) 



com’); {re-execute plot.com) 
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Program Roots; 

label 1; { label decler at ion for goto statement) 

{$! Typedef .sys) 

TYPE 

ary4 = array 
ary 6 = array 

var 



[1 . . 1 1 ] 
[ 1 . . 21 ] 



of real ; 
of integer ; 



posCounter , i , j , j j ,kk1 ,m, 1 . 
cnpoles , sizezeros , one : integer; 
dencoef f ,numcoef f , cdencoeff : ary4 ; 
r ealpartzero , imagpar tzero , 
crealpartpole , cimagpartpole : ary3s; 

PSI : aryls; 

plotfile :file; 
gain :real; 

($I Polynom. inc > , 

{$I Rootf ind . inc ) 

BEGIN 

ClrScr ; one : = 1 ; 

Character i St ic_equation( A1 A , size ,Dencoef f ) ; 
for i;=1 to size, do 



PSlti,size]:=B[i,l] ; 

do 

to size do 



for jj:=2 to size 
begin 
for 



i : = 1 
begin 

1 :=size-jj+1 ; 
klc 1 * = n +1 * 

PSI ti > 0 J : =Dencoeff [kkl ] * B[i,l]; 

for,l:=i to size do 

®fsi[i, j] :=PSI[i,j] +AlA[i,l] * PSI[l,kkl]; 
end; 
end ; 
end; 

for i;=1 to size do 
begin , ^ 

Numcoeff [i J ; =0 . 0 ; 
for J:=1 to size do 
begin 

Numcoeff [i] : =Numcoeff[i] + PSI[j,i] * 
end; 
end; 

for i:=1 to size do 
begin 

m ! =siZG+1 ”i * 

if numcoeff tm] <> 0.0 then goto 1; 
end ; 

1 : sizezeros : =m-1 ; 

root_f inder (sizezeros , Numcoeff , r ealpartzero , 
imagpartzero , one ) ; 

ain:= Numcoeff [sizezeros+1 ] ; 
or i;=1 to sizezeros+1 do 

Numcoeff [i] : = Numcoeff [i] /gain; 
end; 

for i:=1 to maxorder do 

CDenCoeff[ii := 0.0; (initialize) 
for i:=1 to SizeZeros + 1 do, , 

CDenCoeff [i J := Numcoeff [i] * gain; 
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for i:=1 to Size + 1 do 

CDenCoeffLiJ := CDenCoeff[i] + Dencoef f [i] ; 



if Size 
{ NPoles 



> SizeZeros then CNPOLES : =Size 
should always be greater) 

else CNPOLES : =SizeZeros ; 



^ (compute new denominator roots) 

root_f inder ( Cnpoles , CDenCoef f , CReal Part Pole , 
CImagPartPole , one ) ; 



ClrScr 
wr 



rScr ; textcolor ( lightblue ) ; { on-screen 
iteln( ' P^ant Characteristic 

Roots 



* » * > ) 



writeln; writelnrhi 
writelnC ’ ROOTS OF T 



hvideo : 

E NUMERATOR: ») 



titles ) 
Equation 



For I:=1 to SizeZeros do (position for output) 
begin 

PosCounter : = (I mod 2 ) ; 

If PosCounter = 1 then writeln; 

LowVi 
write 

write 
end; 



deo: (write zeros } 

( ’s[’ ,1, ’ J = ’ ,Real 
’ +j ’ , ImagPart^er 



( 






, Real Far tZerofl] : 1 0 
- o[IJ:10:M; 



writeln: writeln: HighVideo: 
writeln( ’ROOTS OF THE DENOMINATOR: ’); 



For I:=1 to CNPoles do 

begin (compute on-screen position) 

PosCounter : = (I mod 2 ) ; 

If PosCounter = 1 then writeln; 



LowVi deo ; 
wr ite ( 



write( ’ 
end; ( for ) 



s[’ 

+d 



^ I C Jma. 



(output roots) 

’ , CReal PARTPOLET iJ: 10: 
gPART POLEl I J : 1 0 : 3 ) ; 



HighVideo ;gotoxy( 1 ,24) ; _ _ _ 

write(’ Press any key to continue or [Shift] [PrtSc] 
for hardcopy.’); 



(check keyboard buffer for value change. 

If number changes by 1 or 2 indicates that shift key 
depressed. If so, then remove "Press any key..." 
prompt from screen so it won’t print to printer) 



keyold := mem[0000 : 1 047] ; not_erased := true; 
Repeat 

key := mem[0000 : 1 047] ; 

. ^ , ^ ^ Q 



if ((key=keyold 



(not erased) then 
begin 

GotoXY( 1 ,24) ; write( 
not erased := false; 
end; 

Until KeyPressed; 



or (key=keyold + 2)) and 



; 80 ) ; 



Assign (plot file , ’ Plot . com ’ ) ; 
Execute (plotfile ) ; 



end. (Program Roots) 
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Program Matrix_Mcinipulation( input , output ) ; 

($I Typedef .SYS) 

($1 Box. INC) 

($I Ut-mod01 . INC) 

var 

help3f ile , inversefile , determinfile , 

Polynomf ile , eigenfile :file; 

Procedure Matr ix_mathemat ics_menu ; 

var 

i,Tab : Integer: 

Okchoice : Set Of* Char ; 

Finished : Boolean; 

Procedure Menultem( Pick:Char ;Description :Str80; 

Color : Integer ) ; 



Begin 

TextColor ( color } : ^ 

Write (’ ^:tab,M’); TextColor (white ) : Write(pick); 

TextColor ( color ); Writeln(’)’ , Be script ion) ; 

End; 



Begin 



6lrScr : TextColor ( white ) ; Finished ;= Fal 
GoToXYC 19,4); 

WriteP*^*^* MATRIX MATHEMATICS MENU 



se ; 



matrix mathematics menu ^^^^»>); 

for i := 1 to 4 do writeln( ’ ’ ); 

Tab := 16; 

MenuItemC’D’ , ’Determinant of A Matrix’ .yellow) ; 
Menultem( ’ C ’,’ Character 1st ic ^polynomial of A 

>y©liow) * 

Menultem( ’ I ’,’ Inverse ’ of A Ma"f;rix’ yellow); 
MenuItemC ’ E ’,’ Eigenvalues of A Matr ix’ .yellow) ; 
Wr iteln ; 

Wr iteln ; 

Menultem( ’H ’ , ’Help ’ , llghtmagenta) ; 

Menultemp Q ’ , ’ Quit to *' ' 

i); 



een 



XJL U lXiQ.pt> W XX UQ / f 

VS Main menu ’, llghtmagenta) ; 



TextColor (gre 

Tex-f;Color?white ) ; GoToXY (40,22); 

Repeat 

Option; 

until ch in [ ’ D ’ , ’ C ’ . ’ I ’ , ’ E ’ , ’ H ’ , ’ Q ’ ] ; 
case ch of 
’D ’ : begin 

Assign( determinf ile , ’ determin . chn ’ ) ; 
Chain( determinf ile ) ; 
end: 

’ I ’ : begin 

Assigni inversef ile , ’ inverse . chn ’ ) ; 
Chain( inversef ile ) ; 
end: 

’ C ’ : begin 

Assign( Polynomf lie , ’ Polynom. chn ’ ) ; 
Chain(^Polynomf ile ) ; 
end; 

’E’ : begin 

Assign( eigenf ile , ’ eigen . chn ’ ) ; 
ChainCeigenf ile ) ; 
end; 

’ H ’ : begin 

Chlin?^iel^?Fiie1- ^ ’ 

end ; 
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» Q ’ : begin 

Assign( SvsFile , ’ svs . COM ’ ) ; 
ExecuteC SvsFile ) ; 

, end; 

End; (case) 

End; 

Begin {main program) 

Clr scr ; 

Repeat 

begin 

Mat r ix_Mat hematic s_Menu ; 
end; 

Until Finished = True; 

End. {main program) 
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Program Matrix Determinant ( input , output ) : 
label 1 0 ; label~20 ; label 30; label 40; label 50; 
label 70;label 80; 



Page 

label 60 ; 



1 



$I 

|l 

$I 



var 



Typedef .SYS 
Ut-mod01 . IN 
Boxuser . INC 



matrixfile :file; 
list :text; 

counter ^ 1 , 3 ,ii,k,m,n, even 



_ tinteger; 

temp value , value , det , det_correction , 
determinant old, determinant :real; 

A1 “ : aryls; 

Begin 

counter 
Boxuser , 
for i:=1 
Begin 

for 3 : = 1 

Begin ^ ^ 

„ A1 [i,3] := A1A[i,3] ; 

End; 

End; 

for i : = 1 
Begin . 
k I — X * 

30:ii‘ A1[k,i] <> 0.0 then goto 10; 



: = 0 ; 

to size do 
to size do 
1 := A1A[i,; 

to size do 



k; =k+1 ; 
if (k-size) 
goto 40; 

|0:^f H-k) 

for 

Begin _ _ 

temp value : =A1 [i ,mj ; 



<= 0.0 then goto 30; 
0.0 then goto 40; 



0:if (, 1 -k} > 0.0 then goto 
f (i-kl = 0.0 then goto 70; 
or m:=1 to size do 
gin _ 

temp value : =A1 [i , 

A1 Li,mJ:=A1 [k,m] ; 

A1 [k , m] : =temp_val 



-ue ; 

End; 

counter := counter+1 ; 

70:ii:=i+1 ; 

if ii > size then goto 20; 
for m:=ii to size do 
Begin 

if Al[m.il = 0.0 then goto 80; 
value:=Al [m,ij / AlLi,lJ; 
for n:= i to size do 
Begin _ _ _ _ _ 

A1[m,nJ:= A1 [m,n] - A1[i,n] * value; 
End; 

80 :End; 

20 :End; 
det : = 1 . 0 ; 

for i:=1 to size do 
Begin 

det:=det * A1[i,i]; 

End; 

det correction:= exp( counter * 
determinant old:=det correction 
even : = counter mod 2T 
if even <> 0 then goto 60: 
determinant := determinant old; 
goto 50; 

60 : determinant : =-determinant old; 
goto 50; 

40 : determinant : =0 . 0 ; 

50 :clrscr ; wr ite In ;TextCol or (yellow) ; 
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Writeln(’The given matrix is 
to size do 
to size do 



); writeln; 



for i : = 1 
Begin 

f or . 3 : = 1 
Begin 

write( ’ ’ , A1 A[i , 3 ] : 1 1 ) ; 

End;Writeln; ^ 

End; Writeln; Writeln; 

Write(’The determinant value is :’)» 

Text Co lor (white ) ; writ eln( determinant ) ; 

keyread(key ,keyold,not_erased) ; 

textcolor (yellow) ; 



R © ID © cL^ 

box pr inter ;gotoxy( 58 , 1 5 ) ; 
writelnf ’ "DETERMIN . RES" ’ ) ; 



gotoxy 



read 






:write( ’ 

) ; 



); 



X w GL\^ \ TV W \mL f X.X J f 

if (ch=’F’) or (ch=’f’) or (ch=’P’) or (ch=’p’) 
then begin^ 

if (ch=’F’) or (ch=’f’) then 
begin 

fOtox^^45 ,15); textcolor ( red^ ; 



write( ’PRINTI 

Assign( list , ’ Determin . RES ’ ) ; 

Rewritet lis-f;) ; 

if {ch=’P’) or (ch=’p’) then 
begin 

otox^^ 45 , 13 ); textcolor ( red) ; 



given matrix is:’); 



write( ’ PRINTI- 

As sign (list, ’LST : ’ ) ; 

Rewritet list ) ; 
end; 

wr iteln( list , ’ 

Result ’ ) ; 

wr itelni list ) : 
write In (list.. The 
wr itelni list ) ; 
for i:=1 to size do 
begin 

for d:=1 to size do 
begin ^ ^ ^ , 

write(list,’ ’ , A1 A [i , 3 ] : 1 1 ) ; 
end;writeln( list ) ; 
end;wr iteln( list );writeln( list); 
Write! list ,’ The determinant value 
writeln! list. determinant ) ; ^ 
writeln( list); writeln! list); 
close! list ) ; 
end; 

Until ch in C’Q’,’q’]; 

Assign(matrixf ile , ’matrix.com’ ); 
Execute(matrixf ile ) ; 



DETERMINANT 



is 



); 



End . 
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Program Matrix inverse( input , output ) ; 
($I Typedef . SYST 
($I Ut-mod01 . INC) 

($I Boxuser.inc) 
type 



var 



ary4 

ary6 



= array 
= ar r ay 



[1 . . 1 1 ] 

[ 1 . . 21 ] 



of real; 
of integer; 



matrixfile :file; 
list ;text; 

S olynom coeff ;ary4; 

inverse : aryls; 

singular , step , stepper ; integer; 

($I Polynom.inc) 

($I Inverse. inc) 

begin 

clrscr; Boxuser ; singular :=0; 

characxeristic_equation(A1 A, size , Polynom_coef f ) ; 




for step;=1 to size 
begin 

for stepper :=1 to size do 
begin , ^ ^ , 

write( ’ ’ ,A1 A[step , stepper] ; 1 1 ) 

end; wr iteln ; 
end; wr iteln ; 
if singular=l then 

writeln(’The matrix is singular 

else 



) 



begin 

writeln(’The inverted matrix is :’) ;wr iteln 

for step;=1 to size do 

begin 

for stepper :=1 to size do 

begin , _ _ 

wr ite( ’ ’ , Ainver se Lstep , stepper] : 1 0 , ' 

end ;wr iteln ; 
end; 
end; 

keyread(key ,keyold,not_erased) ; 

Repeat 



); 



Eox printer ;gotoxy( 58, 1 5 ) ; 
writelnC ’ " INVERSE . RES ; 



write( ’ 



textcolor 
gotoxyC 1 , 



H9.1 j); 



1 1 ow ) ; 

' ); 



read(kbd, ch / , , ^ ^ , 

if (ch=’P’} or (ch=’f’) or (ch=’P’) or (ch=’p’) 
then begin ^ 

if fch=’F’) or (ch=’f’) then 
begin , ^ ^ 

got oxy( 45 , 1 5) ; textcolor ( red ) ; 

writeC ’ PRINTING ’); 

As sign (list, ’ Inverse . RES ’ ) ; 

Rewritet list ) ; 
end; 

if (ch=’P’) or (ch=’p’) then 
begin 
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gotoxyf 45,13); text color (red) ; 

writeC ; PRiNTING 

Assign( list , ’ LST : ’ ) ; 

Rewrites list ) ; 
endj 

wr iteln( list , ’ 

, , RESULT ’ ) ; 

wr iteln( list ) ; 

writelnl list ,* The given matrix 
writeln(list); 

step:=1 to size do 



’ ); 

MATRIX 

is:’); 



INVERSE 



to size do 



for 
begin 

for stepper : = 1 

begin , ^ ^ 

write(list, ’ ! ,A1A[step, stepper] : 1 1 ) 

end;writeln(iist); h 

end;writeln( list);writeln(list); 
begin 

writelni list The inverted matrix is:’) 

wr iteln( list ) ; 

for step:=1 to size do 

begin 

for stepper :=1 to size do 
begin 

wr ite( list , ’ ’ , _ 

Ainver se [step , stepper J : 1 0 , ’ 
end; writ eln( list ) ; 
end; 

end: ^ ^ . 

writeln(list); write In (list); 
close( list ) ; 
end; 

Until oh in [’Q’,’q’]; 

Assign( matr ixf i le , ’matrix. com ’ ) ; 

Execute(matr ixf ile ) ; 



); 



end . 
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Program Matrix Manipulation( input , output ) 
($I^Typedef .SYSi ^ 

{$! Ut-mod01 . INC}{ $I Boxuser.inc) 



type 



ary4 

ary5 

aryo 



array . 1 
array . 1 
array [ 1 



1 1 ] 


of 


real ; 


21 1 


of 


real ; 


21] 


of 


integer 



var 

matrixfile :file; 
list :text; 

c2,C1 :ary4; 

Real Part , ImaginaryPart : ary3s ; 
i,j,one : integer ; 

($I Polynom. inc > , 

{$! Rootf ind . inc ) 



begin 

clrscr ; Boxuser ; one : =1 ; ^ ^ 

character i St ic_equat ion ( A1A,size,C1 ); 

root_Finder ( size , Cl , RealPart , ImaginaryPart , one ) ; 
clrscr ; Text col or (yellow) ; 

writelni’The given matrix is : ’ ) ; wr iteln ; 

for i:=1 to size do 

begin 

for j:=l to size do 
begin ^ , 

write( ’ ’ , A1 A[i , J ] : 1 1 ) ; 

end; wr iteln ; 
end; wr iteln ;wr iteln; 

writeln(’The eigenvalues of the matrix are: ’ )i 

writelnC ’ REAL PART 

IMAGINARY PART ’ ) ; 

for i:=1 to size do 
begin 
write 
write 
write 
write 
end; 



iRealPart [i] : 3 ) ; 
ImaginaryPart; [i] :3 j ; 



); 



wr iteln ; 



keyread(key ,keyold,not_erased) ; 



Repeat 



P®2Lt 

box printer :gotoxy( 58 ,15); textcolor (yel 
writelnC ’ "EIGEN. RES" ’ ) ; gotoxyC 1 ,25) ; 



write( ’ 



if (oh 






* , or (ch=’f’) or (ch=*P’) 

then begin ^ 

if Cch=’F’) or (ch=’f») then 
begin , ^ ^ 

gotoxy(45 , 1 5 ) ; textcolor ( red) ; 

writeC ^ PRINTING 

Assign(list, ’Eigen. RES’ ); 
Rewritet lis-t ) ; 

if {ch=’P’) or (ch=’p’) then 
begin 

gotoxy( 45,13); textcolor ( red) ; 

writeC ’ PRINTING. ... 

Assign( list , ’ LST : ’ ) ; 
RewriteClist ) ; 
end; 



or 



); 



); 



low) ; 

); 

( ch= ’p ’ ) 
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wr iteln( list , ’ 

, RESULT 
list ) : 
list , ’ The 
list ) j 
to size do 



wr iteln 
wr iteln 
wr iteln 
for i : = 
begin 
for .0 
begin 



’ ); 

given matrix 



IS 



Page 

EIGENVALUES 

); 



=1 to size do 



wr ite( list , ’ 
end;write}n(li! 



. ; ,A1A[i,j] :1 1 ); 

1 S "b ) * 

end:wf iteln^^list ) ;writeln(list); 
wri-teln(list , ’ 

The eigenvalues of the n 
writelnC list , ’ 

IMAGINARY PART ’ ) ; 
for i:=1 to size do 
begin 

writeC list , ’ _ _ . 

writeClist , Real Par t [ij : 1 0 ) ; 
writeClist,’ +1 _ ;); ^ 

wr iteC list , ImaginaryPart [i J : 1 0 ) ; 
wr iteln( list ) ; 

writeln( list );writeln( list); 
close( list ) ; 

Until’ch in [ ’ F ’ , ’ f ’ , ’ P ’ , ’ p ’ , ’ Q ’ , ’ q ’ ] ; 

Assign( matr ixf ile , ’matrix . com ’ ) ; 
Execute (matrixfile ) ; 



matrix are 



REAL PART 



); 



h 



end . 
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Program Matrix Characteristic Equation( input , output ) ; 
($I Typedef.SYSK - ^ 

($I Ut-mod01 . INC) 

{$! Boxuser . inc } 
type 



ary4 

ary6 



= array 
= array 



[1 . . 1 1 ] 

[1 . . 21 ] 



of real; 
of integer; 



var 



Cl 



matrixfile :file; 
list ;text; 



: ary4 ; 

: integer ; 

: str ing[2J ; 



i , j , vertpos ,horizpos ,poscounter 
exponent 

($I Polynom.inc) 

begin 

clrscr : Boxuser ; 

characteristic equation(A1 A,size,C1 );clrscr; 
Textcolor (yellow) ; writeln; 
writeln(’The given matrix is : ’);writeln; 

for i:=1 to size do 
begin 

for j:=1 to size do 
begin , _ ^ , 

write(* ’,A1A[i,j]:11); 
end; writeln ; 
end;writeln; writeln; 

writeln(’The characteristic polynomial of given 
matrix is : ’ ) ; 

vertpos : =size+7 ; 
for i:=size+1 downto 1 do 
begin 

poscounter : = ( 3 mod 4) +1 ;horizpos : =poscounter * 

if poscounter = 1 then vertpos : =vertpos+2 ; 
if i <> 1 then 
Id 0^i n 

gotoXY([_hor izpos-1 0,vertpos)^write(C1 [i] :7:4); 
msg( ■ * ■ 



14; 



( norizpos-i 0 , vertpos ) ; w 
S +’ ,horizpos , vertpos ) ; 
str(i-1 :2, exponent); 

msg( exponent ,horizpos+1 ,vertpos-1 ); 



end 
else 
begin 



en 

end; 



gotoXY(horizpos-1 0 .vertpos );write(C1 [i] :7:4); 



keyread(key ,keyold,not_erased) ; 
Repeat 

box printer:. 



writelnf’ "P 




otoxy ( 58 . 15); textcolor ( yellow) ; 
YNOM.RES^' ’ ) ; 



wr i t e ( 



); 



X w CL^L V VjL • J 9 

if (ch=’F’) or (ch=’f’) or (ch=’P’) or (ch=’p’) 
then begin 

if (ch=’F’) or (ch=’f’) then 
begin 

gotoxy( 45,15); textcolor ( red) ; 

writeC ’ PRiNTING 

Assign (list, ’ Polynom.RES ’ ) ; 

Rewr ite( list ) ; 
end; 



); 
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’ ) ; 

POLYNOMIAL 

); 



to size do 
to size do 



if (ch=’P’) or (ch=’p’) then 
begin 

gotoxy( 45 , 13 ) ; text col or ( red) ; 

writeC ’ PRINTING 

Assign(iist, ’LST: ’ ); 

Rewritet list ) ; 
end: 

wr iteln( list , ’ 

RESULT ’ ) ; 

writeln(list) ; 

wr itelni 1 ist . ’ The given matrix is 
wr itelni list ) 
for i : = I t 
begin 

for 3 : = 1 

begin ^ ^ ^ , 

writedist,; ’ ,A 1 A[i ,3 ] : 1 1 ) ; 
endjwr iteln( list ) ; 
end;wr iteln( list);writeln(list); 
write( list The characteristic polynomial of 
given matrix is:’); . ^ 

wr itelni list , ’ ( Descending power of S ) ’); 

writeln( list ) ; 
for i:=size+1 downto 1 do 
begin 

write( list ,C1 [i] :7*d» ’ ’)» 

str ( i-1 : 2 , exponent ) ; 
end : 

writeln(list); write In (list); 
close( list ) ; 
end; _ _ 

Until ch in [’Q’,’q’]; 

Assign(matr ixf ile , ’matrix.com’ ) ; 

Execute(matrixf ile ) ; 



end . 
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Program Control ( input , output ) ; 

<$I,Typedef .SYS) 

{Program type and variable definitions 
($I Ut-mod01 .iNC > (input utility programs 



type 

var 



Aryl = Array [ 1 .. 1 1 , 1 .. 1 1 ] of real; 

Stepping, Steps , step , steppings , 

Tempi , Temp 2 , valuel , value 2 ,m,mplus , 
last rank,k, valuelplus .kplus 
new matrixl ,new matrix2 



big matrix 
list 



; integer ; 

: aryl s ; 
:ary1 ; 

: text ; 

(finds controllability result) 



t l control . inc ) 

I Boxuser . inc ) 

5 egin 

valuel := 1; value1plus:= valuel; 

Clr Scr ; Text col or ( yellow ) ; 

Writeln(’The A matrix is : ’)jwriteln; 

For steps := 1 to size do 

Begin 

For stepping := 1 to size do 
Begin ^ ^ ^ ^ 

Write ( ’ ’ ,AlA[steps, stepping] : 1 1 ) ; 

End; Writeln; 

End; write In ; writ eln ; 

Writeln(’The B Matrix is 
For steps ;= 1 to SIZE do 
Begin 

for stepping; =1 to ni do 
begin ^ _ 

Write( ’ B [steps , stepping] 
end; writeln ; 

End; writeln : writeln ; 

For steps :=1 to size do 
Begin 

For stepping :=1 to ni do 
Begin 

new_matrix2 [steps , stepping] 



’ ) ;wr iteln; 



11 ); 



new_matrix1 [steps , stepping] 



= B[steps , stepping] ; 
= B[steps , stepping] ; 



End; 

End; mplus:= size- 1 ; value 2 ; =value 1 ; 

For steps ;= 1 to mplus do 
Begin 

matrix multipllcation(A 1 A, new matrixl . 

new matrixl . size , size .value2 ) ; 
step := valuel +1; valuel ;= valuel + value2; 
For k ;= 1 to size do 
Begin 

For stepping ;= step to valuel do 
Begin 

steppings := stepplng-step+1 ; 
new matrix2[k, stepping] := 

new matrixl [k , steppings] ; 
End; “ 

End; 

End; 

For steps :=1 to size do 
Begin 
For 



stepping; =1 to valuel do 

g_m9.tr ix[ steps , stepping] 

] i 



Begin 
big ma 

new matr 1X2 [steps , stepping 
“End; 
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: NOT_ERASED ) ; 



End; mplus : = 1 ; 

For steps ;=1 to size do 
Begin 

Matrix_reduction(big matrix , mplus , size , last rank); 
mplus := last rank +1; — 

End; ~ 

If* last rank = size then 

WritelnC’The system is controllable. ’);Writeln; 

If last rank < size then 

Writeln“( ’ The system is uncontrollable. ’);Writeln; 
value 1 : =value1plus; 

KEYREAD ( KEY , KEYOLD , : 
repeat 

Read(kbd, ch ) ; 

If (ch=’F’) or (ch=’f’) or (ch=’p’) 
then Begin 

If (cn=’F’) or (ch=’f’) then 
Begin 

gotoxyC 45,15 ); Text col or (red) ; 

writelnC ’ PRINTING 

AssignC list , ’ control . RES ’ ) ; 

Rewr ite( list ) ; 

End: 

If {ch=’P’) or (ch=’p*) then 
Begin 

gotoxyC 45 ,1 5 ); Text col or ( red) ; 

writelnC ’ PRINTING 

Assign(list, ’LST: ’ ); 

Rewr ite( list ) ; 

End; 

Wr iteln( list , ’ 

RESULT ’ ) ; 

Wr iteln( list ) j 



’ ); 

or (ch= 



P’ ) 



); 



); 



CONTROLLABILITY 



Writelnflist. ' The 
wri teln( list ) ; 

For steps : = 1 to 
Begin 

or stepping := 1 



Plant matrix 
size do 

to size do 



A is 



); 



, r -IS 

Write(list^ ’ ’ ,B[steps, stepping] : 1 1 ) 



Begin 

Write(list, ’ ’ ,A1A[steps, stepping] : 1 1 ) ; 

End;Writeln(iist ) ; , , 

End: Writ eln( list ) ;Writeln(list); 

WritelnC list The input Matrix B is : ’); 

wr iteln( list ) ; 

For steps : = 1 to SIZE do 
Begin 

For stepping := 1 to ni do 
Begin 

Wr ite( list , ’ 

End;writeln( list ) ; 

End; WritelnC list) ; Write In ( list ) ; 

If last rank = size then 

WritelnC list The system is controllable. ’); 
Wr iteln ; 

If last rank < size then 

Writeln~( list The system is uncontrollable. 
Writeln: 

WritelnC list ); WritelnC list ) ; 

CloseC list ) ; 

End; _ _ 

Until Ch in C ’Q’ , ’q’ll s 
AssignC SvsFile , ’ svs . COM ’ ) ; 

ExecuteC SvsFile ) ; 

End. 



); 
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Program Observabi 1 i ty ( input , output ) ; 
label Ijlabel 2;label 3;label 4; 

<$I Typedef . SYS) {program type and definitions) 
type 

Aryl = Ar ray [ 1 . . 1 1 , 1 . . 1 1 ] of real; 



var 



$I 

SI 

$I 



stepping , steps , step , temp 1 , temp2 , 

add no , value2 , last rank : integer; 

matrix ■ 

matrix matrix 
list 

Boxuser . Inc ) 

Control . inc ) 

Ut-mod01 .inc) 



aryl s ; 
aryl ; 
text ; 



(utility input routine) 



Begin 

ClrScr ; Text color (yellow) ; 
Writeln(’The Plant matrix A i 
:= 1 to size do 



s : ’ ) ; wri teln ; 



For steps 
Begin 

For stepping := 1 to size do 

Begin , _ _ ^ 

Write( ’ ’ ,A1A[steps, stepping] : 1 1 ) ; 

End;Writeln; 

End; Wri teln; wri teln; 

Writeln(’The output Matrix C is : *);>^riteln; 

For stepping := 1 to no do 
begin 

for steps :=1 to size do 
begin 

Write( ’ C [stepping , steps] : 1 1 ) ; 

end; wr iteln ; 
end: 

wr iteln ; wri teln ; 
add_no;=no; value2 ;=1; 

For steps := 1 to no do 
Begin 

For stepping: =1 to size do 
Begin _ _ 

matr ix_matr ix[steps ^steppin^j := 

matrix [steps , stepping^ ;2 steps , stepping] ; 
End; 

End; 

4:matrix reduct ion( matr ix matrix,? 

size, last rank 
(calculate the rank of i 
If last rank < size then GoTo 2; 

Writeln(’The system is observable with index -, 

value2,’ );Writeln; 

KEYREAD(KEY,KEYOLD,NOT ERASED); 

GoTo 1 ; 

2: If value2 < size then GoTo 3; 

Wr iteln (’The system is unobservable. ’); 

KEYREAD ( KEY , KEYOLD , NOT ERASED ) ; 

Writeln:GoTo 1; 

3 t VS.l'LlGiC t = V3 l1\JG2+1 * 

matrix mui tipi i cat ion (matrix, A 1 A, mat rix.no , 

SXZ0 siz©)* 

(multiply the A and C matrices ) 

For steps := 1 to no do 
Begin 

step:= last_rank + steps; 
for stepping: =1 to size do 
Begin 



idd no , 

k • 

:Ae matrix) 
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matr ix_matrix[step ^stepping] := 

_ ^ matrix[steps , stepping] ; 

End; 

End; 

add no := last rank + no; 

GoTU 4; - 

1 : repeat 

Box printerigotoxyC^ 5) ; Textcolor (yellow) ; 
writeln( ’ "OBSER.RES" ’ ) ; gotoxy( 1 , 25 ) ; 
wr ite( ’ =. O' ’ ) : 

§otoxy(49, 17) ; 

Read(Kbd, ch ) ; 

If Cchs’P’) or (ch=’p’) or (ch=’F’) or (ch=’f’) 
then Begin 

If (cn=’F’) or (ch=’f’) then 
Begin 

GotoXY( 45,15 );textcolor ( red) ; 

writelnt ’ PRINTING ’); 

AssignClist, ’obser.RES’ ); 

_ . , {print , observability output ) 

Rewr 1 te ( 1 ist ) ; { to file on the current drive) 
End; 

If (ch=’P’) or (ch=’p’) then 
Begin 



gotoxyC 45,13 ); textcolor ( red) ; 
writelni ’f'RIN'tiNG 



); 



AssignClist , ’ LST : ’ ) ; 

{print observability result to the printer) 
Rewr ite{ list ) ; 

End; 

WritelnClist , ’ OBSERVABILITY RESULT 



); 



Writeln( list) ;Writeln( list) ; 

Writelni list The Plant matrix A is :’); 
wr itelnC list ) ; 

For steps ; = 1 to size do 
Begin 

For stepping : = 1 to size do 
Begin 

Write(List, ’ ’ ,A1A[ steps, stepping] : 1 1 ) ; 

End;Writeln(iist); ^ sj 

End ;Writeln( list ) ; Writeln(_list ) ; 

Writeln( list The Output Matrix C is : ’); 

writelnC list ) ; 

For steps := 1 to no do 
begin 

for stepping: =1 to size do 
begin 

Write(list,’ ’ , C [steps , stepping] : 1 1 ) ; 
end ;writeln( list ) ; 
end: Writeln(List 7 ; WritelnC 1 ist ) ; 

If last rank = size then 

WritelnC list ,’ The system is observable with 
, index ’,value2,’ .’); 

If Cvalue2 >= size ) and Clast_rank < size) then 
Writeln (list, ’The system is unobservable. ’); 
writeln(list); WritelnClist); 

Close(list7,* 

End; _ _ 

Until Ch in C’Q’ , ’q’ ] ; 



AssignCSvsFile, ’svs.COM’ ); 

{re-execute SVS main program) 
ExecuteC SvsFile ) ; 

End. 
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Program Luenberg observer ( input , output ) ; 
(SI^Typedef .SYS>- 

($i Ut-mod01 . INC} 

{$! Ut-mod02 . inc > 

label 940;label 960;label 22;label 28 ; label 950; 



Type 



Var 



aryl 
ary2 
ary3 
ary 4 
ary 6 



array 

array 

array 

array 

array 



1 . . 1 0] of real ; 

1 . .30. 1 . .30] of real ; 
1 . . 30. of real ; 

. 1 . . 1 1 . of real ; 

1 . . 1 1 ] of integer ; 



rrl . rim, ooo , oooo 

i,NRl,l ,mm,nrm,nrp , ixx, jxx, s , one 
nrn , nr , r , rp , k , t , o , umr an , i i , j j , 
code , ver tpos , hor izpos , poscounter 
h , a,1 , e 

coeff , coef , desired feedback 

f “ 

beta,x 

phi , u 

list 

strg , exponent 
change 

specification 

Realpartvalue , imaginarypartvalue 
temp , input type 



real ; 



: integer ; 
aryl ; 

: ary4 ; 

: aryl s 
;ary3 ; 

: ary2 ; 
;text; _ _ 

; str ing[2] ; 
iboolean; _ 

: str ing[5] ; 
; ary3s ; 

: char ; 



{$! Ut-mod03 . inc ) 

{$! Luenberg . inc ) 

{$! Pole. inc) 

\CI Rootf ind . inc ) 
($I Boxuser.inc) 



Begin 

ClrScr ; luen ; =true ; one : = 1 ; 
for i:=1 to size do 



do 



' ); 



)rl: = 1 to 1_d 
E[i1;=Bri,l]; 
gotoxyC 1,22); 

invvideoC ’ Press <ESC> to change it!, 
gotoxyC 1 , 23 ) ; 

invvideoC ’ Then type your input with <ENTER> key 
gotoxyC 1,1); Textcolor C lightblue ) ; 
writelnC ’ *** Luenberger Observer Design 

Parameters 
TextColor Cyellow) ; 

writelnC ’ = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = ’ ) ; 
Textcolor Cyellow) ; 

msgC ’ Input degree of observer C 10 max)’ ,1,5); 
Repeat 

inpytC ’N’ , ’ ’ , 40 , 3 , true , F1 ,F10); 



val ( answer , r , cod 



rr ^ 1 *1 p 

10) and 
<= 10) anc 






< 1 ) 

> 0 ); 



then beepC 900 , 350 ) ; 



if Cr > 
until Cr 
rp ; =r + 1 ; 

msgC ’ Input the Desired Feedback Coefficients in 
, Factored <F> Form ’,1.7); 
msgC ’ or Coefficient <C> Form ’,1,8); 

repeat , 



Input C ’A’ , ’C’ ,64, 7, 2, true, FI ,F10); 
temp : =copy C auiswer ,1,1); 
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C;F’ |C’]) then beep(900 ,350) ; 
, ’ C ’ J ; 



if not (temp^in 
until temp in [’F 
inputtype : =temp ; 

Case inputtype of 
’ F ’ :begin 

if change then 
for i:=1 to size do 
begin 

str f realpartvalue [i] :8;2^Filvar [2*i+9] ) ; 
str ( imaginarypartvalue [i1 : 8 : 2 , ^ 

Filvar[2»i+10] ); 

end; 

input_Factored( ’POLES’ , size , Real Partvalue , 

ImaginaryPartvalue ) ; 

end: 

’ C ’ rbegin 

if change then 
for i;= 1 to size-1 do 
begin, 

str(Desired f eedback[size-1 ri] :8 : 2 , 

“ Filvar [size+22-i] ) ; 

end; 

input coef (’ POLES ’, size-1 , Desired feedback) 
end; 

end ; Clr Scr ; Wr iteln ; gotoxy ( 1 , 22 ) ; 

invvideo( ’ Press <ESC> to change it!, ’ 

gotoxy ( 1 , 23 ) ; 
invvideo( ’ Then t, 
sotoxvv 1 1 / * Tsx 

writeln( ^ Luenberger Observer Design 

Input Parameters 
textcolor ( yel low ) ; 

writelnC ’=========================================’); 

"t ©X”t c o 1 o r* ( © ©n ) ^ 

msgC ’ Input observer characteristic polynomial 

Factored (F) Form ’,1,5 
or Coefficient (C) Form ’ 



); 



e your input with 
olor ( lightblue ) ; 



< ENTER > key ’); 



msg( ’ 

r ©ID © 3.~t« 

input ( ’A’ , ’C’ ,65, 5, 2, true, FI ,F10); 
temp:= copy ( answer ,1 , 1 
" f not (temp in [ ’ 

il temp.in r.’F’ , ’C’ J ; 



f ' f 



6 ); 



1 

until temp 
if temp = ’C 



hen 



^]) then beep( 900 , 350 ) ; 
goto 950; 



if change then 
for i ; = 1 to r do 

begin , ^ ^ ^ 

str ( realpartvalue[i] ;8i2.Filvar [2’»i+9] ) ; 
str ( imaginarypartvalue [ij : 8 : 2 , Filvar [2*i+10]); 

end; 

input Factor ed( ’POLES’ , r , Real Partvalue , 

Imaginarypartvalue ) ; 

Polynomial of roots( r ,Realpartvalue , 

~ ~ Imaginarypartvalue , coef ) ; GoTo 960; 



950: if change then 
for i : = 1 to r do 



^^^^str ( coef [r-i] :8; 2 , Filvar [r+22-i] ) ; 
end; 



input_coef f ( ’ POLES ’ , r , coef ) ; 

root f inder ( r , coef , realpartvalue , ^ 

“ imaginarypartvalue , one ) ; 
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: =0 . 0 ; 
r do 



r do 
r do 



ClrScr ; Writ eln; Writ eln; 

Writeln(’The observer eigenvalues are : 

WritelnC ’ REAL PART IMAGINARY pArT ’ ) ; 

for i : = 1 tor do 

Begin , _ _ 

Write( ’ ’ , realpactvalue [i] : 1 0 , ’ + ’ 

imaginarypartvalue [i ] : 1 0 , 3 ’); Writeln; 

End: 

Wi*ix©ln *Wri"t©Xri * 

Writeln{’The Observer characteristic polynomial 
coefficints in ascending powers of S’); 

Writeln; 

for i :=1 to rp do 
Begin , ^ ^ 

Write( coef [i] : 1 0 , ’ ’); 

End; Writeln ; Writeln : 

960 :KEYREAD(KEY,KEYOLD,NOT ERASED ); ClrScr ; 
for i:= 1 to 30 do ~ 

Begin 

for 1 : = 1 to 30 do 
Begin ^ 

^ |HI[i,l] 

End; 

End; 

for i : = 1 to 

0 . 0 ; 

§?fi:=1.0: 

for 1 := 1 to 
Begin 

for 1 : = 1 to 
Begin 

Fti,l] := 0.0; 

Eijd; 

F 
i 
F 

28 : End ; 

for i:=1 to 20 do 
begin 

BETiti?:26.0; 

End; 

nr := size * r; 
for i : = 1 to nr do 
Begin _ _ 

BETA[i] := 0.0; 

End; 

nrp : = nr+1 ; 
nrn := nr +size; 
for i : = 1 to size do 
Begin ^ _ _ _ 

BETA[nr+i] := desired feedback[ij ; 

End; 

for i : = 1 to r do 
Begin 

ii:= size*(i-1); 
for 1 : = 1 to r do 
Begin 

: = size * ( 1-1 ) 
for k:= 1 to size 

®®?HI[ii+k,33+k] 

End; 

End; 

End; 

for ii:= 1 to r do 



^ i,l] := -COEF[rp-il; 
.:r 1 = r then goto 28; 
^ i ,i+1 J : = 1 .0; 



do 

= -F[i,l] ; 
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* size; 
to size 



1 



do 

to size do 
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for 1:= 1 
Begin 

for k : = 

Begin ^ 

PHI[i+l,i+k] ;= PHI [i+1 , i+k] + A1A[k,l]; 
End; 

End; 

End; 

for i i : = 1 to r do 
Begin , , . . . v 

1 := size * (ii-l); 
for 1:= 1 to size do 
Begin ^ 

PHI [nr+1 , i+1] := H[ii] ; 

End; 

End; 

for i:= 1 to size do 
Begin 

for 1:= 1 to no do 
Begin^ 

Phi [ nr+i , nr+1] 

End; 

End ; 

nrm := nr +no ; 
for ii : = 1 to r do 
Begin , .... 

1 : = size * ( 11-1 ) : 
t;= nrm + no *(ii-1); 
for k:= 1 to no do 
Begin 

for s:= 1 to size do 
Begin ^ 

PHI [i+s , t+k] 

End; 

End; 

End; 

mm:= r* (size+no) + no; 



= C[1 , i] ; 



= -C[k,s]; 



Boxuser ; ^ 

Text Co lor (white ) ; 

linear_equation(PHI,NRN,MM,BETA,X,K,U) ; 

for ii:= 1 to r do 
Begin 

i : = size * ( ii-1 ) ; 

AJ[iij := 0.0; 
f or o := 1 to size do 

^^J^ii] : =AJ [ii] +E[o] *X[i+o]; 

End; 

End; 

nr 1 : = nr + 1 ; ^ ^ 

ClrScr ;Wr iteln ; Wr iteln ; text color (yellow) ; 
Writeln(’The F Matrix is : ’);Writeln; 

for i : = 1 tor do 
Begin 

for o:=1 to r do 
Begin 

write( ’ ’ ,F [i , o] : 1 1 ) ; 

End ; Wr iteln ; 

End; 

KEYREAD(KEY,KEYOLD,NOT ERASED); ClrScr : Wr it eln : 
textcolor (yellow) ;Writeln( ’ The G1 Matrix is : M; 

Wr iteln ; 

for ii:= 1 to r do 
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Begin 

i:= no + nrm; 

for o:= 1 to no do 
Begin 

Write( ’ ’ , X [o+i] : 1 1 ) ; 

End;Writeln; 

Endtwr iteln ; wr iteln; 

Writeln (’The G2 Matrix is 
for i := 1 to r do 
Begin ^ ^ ^ , 

writeln( ’ ’,AJ[i]:11); 

End; 

KEYREAD ( KEY , KEYOLD , NOTERASED ) : 

ClrScr ; text col or (yellow) ; Write In; 

Writeln (’The output feedback coefficients are 
( Ascending powers of S ) :’);Writeln; 
for o:= nrl to nrm do 
Begin 

wr ite( ’ ’ , X [o] : 1 1 ) ; 

End; Writeln ; writeln ; 

Wriieln(’The compensator feedback coefficients are 
(ascending powers of S ) :’); 

writeln;for i ; = 1 to r do 
Begin 
Write( ’ 

KEY6eAD ( KEY , KEYOLD , NOT_ERASED ) ; 

Repeat 

Box printer ; text color (yellow) ;gotoxy ( 58,15); 
write( ’ "LUEl^BERG.RES" ’); 

gotoxyC 1 , 25 ) ; wr ite( ’ ’); 

’=’f”) or (ch=’p’) or (ch=’F’) or (ch=’f’) 



,H[i]:11 ); 



if (ch 
Begin 

if (ch= 
Begin 



F’ ) or (ch= ’f ’ ) then 



gotox.Y( 45,15); Text col or ( red) ; 

write( ’ PRINTING 

Assign( list . ’ luenberg . RES ’ ) ; 



); 



Rewritet lls't ) ; 

End; 

if (ch=’P’) or (ch=’p’) then 
Begin 

gotoxy( 45,13); text col or ( red ) ; 

write( ’ PRINTING ’) 

Assign( list , ’ LST : ’ ) ; 

Rewritetlisi) ; 

End; 

Writeln( list , ’LUENBERGER OBSERVER RESULT 
For i:=j to 2 do wr iteln( list ) ; 

Writeln( list ,’ The plant matrix A is : ’) 

Wr itelni list ) j 

for i:=1 to size do 

Begin 

for 1:=1 to size do 
Begin ^ ^ , 

Writedist ’ ’,A1A[i,l]:11); 

End; Wr iteln( list ) ; 

End;Wr iteln( list ) ; . 

Writeln( list ,’ The input matrix B is :’) 

Wr i telni list ) ; 

for i:=1 to size do 

Begin 

for 1 : = 1 to 1 do 
Begin 

Write(list,’ ’ , B[i , 1] : 1 1 ) ; 



); 



then 



1 1 0 
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End;Writeln( list ) ; 

End ;Writeln( list ) ; 

Writeln( list The output Matrix C 
Wr itelnC list ) ; 

i:= 1 to no do 



IS 



); 



to size do 



for 
Begin 

for 1:=1 
Begin 

WriteClist ’ ’ , C[i , 1] : 1 1 ) ; 

End; wr iteln( list ) ; 

End ;Writeln( list ) ; 

Wr iteln( 1 ist , ’ The desired feedback coefficients 
„ are : ’ ) ; Writeln( list ) ; 

for i:=1 to size do 
Begin 

for 1 : = 1 to 1 do 
Begin 

write(list,’ desired feedback[i]:11); 

End; Wr iteln^iist ) ; “ 

End ;Writeln( list ) ; 

Writelni list The observer eigenvalues are : 
WritelnC list , ’REAL PART IMAGINARY PART’); 

for i : = 1 tor do 
Begin 



Wr ite( list , ’ ’ , realpartvalue [i] : 1 0 , ’ 
imaginarypartvalue [i J : 1 0 , ’ j ; 
End;Wr iteln(list ) ; 



^r iteln( list ) ; 



WritelnC list ,’ The Observer characteristic 
polynomial coefficints in ascending powers of S’); 

Wr i teln( list); 
for i : = 1 

Begin , _ _ , 

Wr i te( 1 ist , ’ ^ ’ , coef [i 1 : 1 1 ) : 

End; Wr i-^eln(. list ) ;_writeln( list ) ; 



to rp do 



’ );Writeln(list); 



WritelnC list ,’ The ' F Matrix'is : ^ ; WritelnC list ) ; 

for i:= 1 to r do 
Begin 

for o:=1 to r do 

Begin ^ ^ 

WriteClist,’ SF[i,o]:11); 

End; WritelnC list ) ; 

End; WritelnC list ) ; 

Wr itelnC 1 ist , ’ The G1 Matrix is : 
for ii:= 1 to r do 
Begin 

i:= no *(ii-1 ) + nrm; 
for o:= 1 to no do 
Begin 

WriteClist,’ ’,X[o+i]:11) 

End; WritelnC list ) ; 

End; WritelnC list ) ; 

Writeln (list, ’The G2 Matrix is 

=1 to r do 



’ ); WritelnC list ) ; 



for i 

Begin ^ , 

WritelnC list ^ ’,AJ[i]:11); 

End; WritelnC list); 

Writeln Clist,’The output feedback coefficients 

are :’); WritelnC list ) ; 
for o:= nrl to nrm do 
Begin , _ _ ^ 

WriteClist,’ ’,X[o]:11); 

End; WritelnC list) ; wr itelnC list 1; 



WritelnC list 
coe 

for i : = 1 tor 
Begin 



,’ The ' compensator feedbaclj 
fficients are :’);Wri 



do 



telnC 1 ist ) ; 
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Write(list,’ ’,H[i];11); 

End; 

for i;=1 to 3 do wr iteln( list ) ; 
Close(list); 



Until oh in E’Q’ 
Assign( SvsFi le , ’ 
Executed SvsFile ) ; 



End; 




e 



End. 
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)rogram optimal , control ( input , output ) ; 
.1 Typedef . sys > 

:i Graphix.sys) 

$I Kernel. sys} 

$I Windows . sys ) 

$I Polygon. hgh) 

$I AxisThgh) 



:i Ut-mod01.INC 
;i Ut-mod02.INC 
:i Ut-mod03.inc 



{ I/O procedures) 



Gr apmenu . inc ) 
" ■ : > 



($I . 

( $I Boxuser . inc 
label 2;label 3; 
const 

xarray : string[14] 

var 

GTNY : array 
GTN ; array 
PSI , PI ,D3,D4 
^D2,^D^ ,GT 
ti , i , 0 ,1 ,nterm,kk , cpde 



y [ . . 

y [ 1 • • 
,FI,q, 



=’TIME INTERVALS 

80] of real; 

10,1.. 80] of real; 
gamma : aryls; 

: ary3s ; 
integer 



y ,ymax ,ymin , r , denominator ,negatif , si .result 

g rapharray .grapharray 1 ;plotarray ; 

umpgraph , qui t :boolean; 

list ;text; 



: real ; 



ans :char; 

{$! Optimal. inc) 

Procedure Pr intGr aphData; 

{this procedure dumps optimal result data to printer) 
Begin 



LeaveGraphic ; Clrscr ; Textcolor (yellow) ; 

Center! PROQAM OUTPUT OPTIONS * * * ^ , ^ , ^ Q ,8Q) ; 

, 1 , 12 ); 



TextColor (green^ , 

msgC ’ Press <P> Print results to the printer 
msg( ’ <F> List results to file name 

"OPTIMAL. RES" on the current dr ive ’ , 1 , 1 3 ) ; 
msg( ’ <Q> Quit ’,1,14); 

repeat 



jpeat^ 

Read(kbd,ch) : 
If (ch = ’FM 



or 

or 



(ch = 
(ch = 



i'] 



F ’ ) or ( ch = 



begin 

if (ch = 
begin 

Assign ( list, * Optimal . RES 
Rewrite! lis-f; ) ; 
end 
else 
begin 

assign( list, ’ LST : ’ ) ; 
rewrite( list ) ; 
end: 

Wr iteln( list ) ; 
write! list ) ; 
wr iteln( list 



or ( ch = 
then 

) then 

’ ); 



P’ ) 



); 



o’ OPTIMAL CONTROL RESULT 
wr iteln( list ) ; . 

write( list ,’ The order of the system is:’); 
write Inf list,size:2); 

wr ite( list ,’ The member of time intervals is; 
write. 



n( list , si : 2 ) ; 
list, ’’The scg.ler 



); 



wr i te . ^ ^ 

writeln(list,r:7:4); 



R is ; ’ ) ; 
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interval is 



); 



wr ite( list The s^ple 
writeln(list,ti : 4 ) ; 
wr iteln( list ) : 
wr iteln( list , The A matrix is : ’ ) ;wr iteln( list ) ; 
1 to size do 



= 1 to size do 



for 1:= 
begin 
for . 3 ; 

begin , ^ ^ , 

write(list,’ ’ , A1 A[1 , 3 ] : 1 1 ) ; 
end;wr iteln( list ) ; 
end: wr iteln( list ) : 

writeln( list The B matrix is : ’ ) ;wr iteln( list ) ; 

for 1:= 1 to size do 

begin 

for 3 : = 1 to ni do 
begin , ^ , 

write(list, ’’,B[1,3]:11); 
end;wr iteln( list ) ; 

end: wr iteln( list ) : , 

writeln( list The Q matrix is : ’ ) ;writeln( list ) ; 
for 1 := 1 to size do 
begin 

for 3:= 1 to size do 

begin , . 

write(list,’ ’ , Q [1 , 3 3 i 1 1 ) J 
end; wr iteln( list ) ; 
end: wr iteln( list ) : 

writeln( list The Pi matrix is : ’ ) Jwr iteln( list ) ; 

for 1:= 1 to size do 

begin 

for 3:= 1 to size do 
begin , _ ^ ^ 

writedist,’ ’ , FI [ 1 , 3 ] : 1 1 ) ; 
end;wr itelnC list ) ; 
end: wr 1 ist ) : 

wr i telni 1 ist . ’ The GAMMA matrixis:’); 

wr iteln( list ) ; 

for 1:= 1 to size do 

begin 

for 3 : = 1 to ni do 

begin ^ ^ ^ 

writeClist , GAMMA[ 1 , j ] : 1 1 ) ; 
end;wr iteln( list ) ; 
end; wr iteln( list ) ; 

writeln(list. ’MINIMIZATION OVER ALL STAGES ’); 

wr iteln( list ) ; 

wr ite( list , ’ N (stages)’); 

for 1:=1 to i do 

write(list,’ GAIN:’,!.’ ’); 



iteln( list ) ; 

ite( list , ’ ’ ) 

r 1 : = 1 to i do 

wr ite( list , ’ 

wr i teln( 1 1st ) ; 

' to ti do 



wr 

wr 

for 



> -j . 2 ’ 

. j • ^ » 

do 



for 3 : = 1 
begin 

wr ite( list , ’ 
for 1 : = 1 to i 
begin 

wr ite( list , ’ ’ 

end;writeln( list ) ; 
end ; 
end; 

until ch in [’F’,’f’,’Q’ 
EnterGraphic ; 
swapscr een ; 



); 



); 



graphArrayl [3 , 1 ] : 1 0 ) ; 
’q’ , ’P’ , ’P’] ; 
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close( list ) ; 
end; 



begin 
cir scr 
writel 

Textcolor* 
writelnf ’ 

'otoXYCi , 25 ) ; InvVideo( ’Press 



; Text Co lor ( lightblue ) ; 

n( ’ *** Design or Optimal Control 

nrocedure ***>); 

(yellow) ; 

deo(’ Press <ESC> to change itT 
hen Type vour input with <ENTER> key’); 
TextColor (green ) ; 



SUM 



Msg( ’ Input number of time intervals for 
, procedure? ( MAX 99 ) ’,1,4); 
input ( ’N’ , ’50’ ,65,4,6, true , FI , f1 0 ) ; 
val ( answer , ti , code ) ; 

Msg(’What IS your sample interval? ’ 
input ( ’N’ , ’0.1 ’ ,35.6, 6, true, FI ,F10); 
val ( answer , si , code ) ; 

Msg ( ’What IS the value of scalar R ? 

Input ( ’N’ , ’0.0’ ,35, 8, 6, true, FI ,F1 0) ; 
val ( answer , r , code) ; 

Msg(’For the following options Which 
function do you want ? ’,1,10); 
TextColor ( lightmagenta) ; 

Msg( ’ COST=terminal+fuel+tran ectory or 

, COST=terminal+trao ectory <0>’,1,11); 
Msg( ’ COST=terminal+fuel 
COST = terminal 



1 , 6 ) ; 

’ ,1 . 8 ); 

cost 



or 

< 1 > ’ , 1 , 12 ); 



’ ’ , 1 , 1 3 7 ; textcolor(yellow) : 

’Where ,1,14); textcolor ( 1 ightmagenta) ; 

’ _ ! » 1 » 1 



Terminal = 1/2 
Trajectory = 1/2 



msg( 
msg 
msg( 
msg( 
msg 
msg 
msg 
msg( 
msg( 
msg 
msg 
msg( 

Input ( ’N’ , ’0’ ,65, 10, 2. true, FI ,F10); 
val ( ajiswer , nterm , code ) ; clrscr ; 
wr iteln( ’ enter the element of Q matrix: 
for i:=1 to size do 
begin 

or j : =1 to size do 



Fuel = 1/2 



X (N) Q X(N) 

-1 T 

X (k) Q X(k) 

N-1 T^ ^ 

U (k) R U 



)! 












readl] 



’ ) ;write( j ) ; 



is 



); writeln; 



«Liid 
write 
write. 

end;writeln; 
end; 

Repeat 

ClrScr ; Writeln: 

Writeln(’The Q Matrix 
for i :=1 to size do 
Begin 

for j := 1 to size do 
Begin , _ ^ ^ 

Write( ’ ’, Q[i,j]:11); 

End; Writeln; 

End; writeln; 

Write (’Do you want to change^ any element 
the Matrix ? ( Y / N ) ’); 

Read(Kbd, Ans ) ; writeln; ^ ^ ^ \ 

(allows user to change entered data > 



of 



); 
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if ( Ans = ’ Y’ ) 

Begin 

write( ’Input the row to 
write(’ Input the column 
readlnC .1 ) ; wr iteln ; 



Program Listing 
or ( Ans =’y’) then 



“eadln( .i ) : wr iteln ; 
rrite{ ^0( ^ ,i ’ , ’ , 

. J J ' =result ; 

C’N’ , ’n’] ; 
-te ) ; 



) ; readln( i ) ; 
: ’ ) ; 

) ; readln( result ) ; 



chan 
to c 



;e : 
icuige 



St 

End; 

Until Ans in 
text CO lor ( whi 
Boxvis©r * 

FI ^d 6 aMMA(SI,FI,GAMMA) : 

( c511 procedure to calculate fi and 

gamma matrices 



TextColgr^gellow) ; clrscr ; 



le FI matrix 
to size do 



is 



); 



size do 
,FI[i,o] : 11); 



wr itelni 
for i : = 1 
begin 

for j : = 1 to 
begin 
write( ’ 
end; wr iteln ; 
end; wr iteln; wr iteln; 
writeln(’The gamma matrix is: ’ 

for i:=1 to size do 
begin 

for j : = 1 to ni do 
begin ^ ^ , 

wr i t e ( ’ ’ , GAMMA [ i , J ] : 1 1 ) ; 

end; wr iteln ; 
end ;writeln;wr iteln; 
keyread(key ,keyold,not_erased) ; 
Boxuser ; 

for i:=1 to size do 
be 



); 



1 ? i ] : =0 . 0 ; 
2 [i1 : =0 . 0 ; 
“ 3:=1 - 

3?i,.i:i 

1 [i . J ] 



D2 

for 

beg 



3 : = 1 to size do 



= 0 . 0 ; 
= 0 . 0 ; ^ 
=Qti .d] ; 



D4 
PI 
end; 
end: 

y : =0 . 0 ;ymax : =0 . 0 ;ymin : =0 . 0 ; 

for kk:=l to TI do 

begin 

denominator : =0 . 0 ; 
for i:=l to size do 
begin 

for d:=1 to size do 

®§1?i] :=D1 [i] + GAMMACd 
end ; 
end; 

for i:=l to size do 
begin 

for d:=1 to size do 

^®§2?i]:=D2[i] + D1 [j] * 
end; 

denominator : ^denominator 
end; 

denominator : =denominator+ 

negatif : =- 1 . 0 : 

for i:=1 to size do 

®§T?i]:= negatif » D2[i] / 



1] » PI Cd ,i] ; 



FiCd ,i] 
+ D1 [i] 
R; 



» GAMMA [ i , 1 ] ; 



denominator ; 
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[I,KK] :=GT[I] ; 



GTN[I ,KKJ : = 
Dill! :=0.0; 
D2 [i J : =0 . 0 ; 



to size do 
to size do 



* PI [1, j] ; 



end; 

for i:=1 
begin 

for . J : = 1 
begin 

PSI[i,j] :=FI[i,j] + GAMMA[i,l] » GT[j]; 
end ; 
end; 

for i:=l to size do 
begin 

for j:=1 to size do 
begin 

for 1;=1 to size do 

: =D3[i , j ] + PSI[l,i] 

end; 
end; 
end ; 

for i:=1 to size do 
begin 

for d*=1 to size do 
begin 

for_l:=l to size do 

®§4?i,d] :=D4[i,d] +D3[i,l] » PSI[l,d]; 
end; 

if^nterm <=^0 then 

PI Ii^.1 J :^D4li.d J t 

3 :D4 [i , d ] *• = 0 • 0 ; 
end; 
end; 

for i:=1 to size do 
begin 

for d*=l to size do 

^^§3?i , d ] : =0 . 0 ; 
end; 
end; 
end; 

for i:=1 to size do 
begin 

for d • = 1 to ti do 
begin ^ ^ _ _ 

STNY[J] :=GTN[I, Jl; ^ ^ 

if abs(y) < T.0e07 then y:=GTNY[dl 

else y : = 1 . 0e07 ; 
if y > Ymax then Ymax:=y; 
if Y < Ymin then Ymin:=Y; 



goto 2; _ 

K ^ 

+ Q[i . a J + 



R 



GT 






oto 3 ; 
* 



Gr aphar r ay C d ,1 
grapharray [i ,2] 

g rapharrayi [d , i J 

; 

Ymax:= 1.2 * Ymax; 

initgraphic ; , ^ 

selectwindow( 1 ) ; 



gotoxy 



3 ; 

y ; 
y ; 



(4,4 


1 ;wr ite< 


’G’ ) ; 


(4,5 


• ;wr itel 


’ A’ ) ; 


(4,6 


1 ;writel 


’ I ’ ) J 


(4,7 


>; write! 


’N’ h 


5,9 


»; write! 


1:2); 



1 17 



FILE: OPTIMAL. PAS 



Page 6 



Prograjn Listing 

drawtext( 250 ,195,1 ,xarray ) : 

Niceaxes( 0 , ti ,ymin ,ymax. ’ ’); 

Selectworld( worldndxglb ) ; ^ 

Selectwindow( windovmdxglb ) : 
Drawpolygon(grapharray ,1 ,-(ti-1 ), 0,0,0); 

repeat until keypressed; 
qui t ; =f alse ; 
r* 013 © 3 ,'t 

Graph menu( ’OPTIMAL CONTROL GAIN 

PLOT * , dumpgr aph , quit ) ; 

if dumpgraph then pr intgraphdata; 
until quitj 
leavegraphic ; 
end ; 

Assign( SvsFile , ’svs.COM’ ); 

Execute( SvsFile ) ; 
end. 
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program pole_placement ( input , output ) ; 

{$I Typedef . sys ) , (comman type & variable definitions) 
{$I Ut-mod01 .inc) 

($I Ut-mod02 . inc > 

label 100;label 110;label 120;label 130;label 140: 
label 150; label 160; label 170; label 180; label 190; 
label 200;label 210; 



type 

ary6 

ary4 



= array 
= array 



[ 1 . . 21 ] 
[1 . .1 1 ] 



of integer ; 
of real; 



var 

AA,PIN,P10 :ary1s; 

CCC :ary4; 

CC,H,DEN,E,HH :ary4; 

realpart , imaginarypart , real_part , 
imaginary_part , realroot . imaginaryr oot , 
realpart 1 , imaginarypart 1 :ary5s; 

nn ,umr an , 1 , j j ,k , 1 , n2 ,n1 ,m ,er s, 
vertpos ,horizpos ,poscounter ,code 
gaini .test 
input2 

exponent , strg 
Change , factored 
specification 
inputf ile 

reduceor der ,msize , one 
list 



: integer ; 

: r eal ; 

: char ; 

: str ing[2] ; 
: boolean : 

: stringCb] ; 
: f lie ; 

: integer ; 

: text ; 



:i Ut-mod03 . inc } 

;i Pole. inc) 

!I Polynom.inc) 

:i Inverse. inc) 

;i Rootf ind. inc ) 

>1 Boxuser.inc) 

begin { open-loop calculations) 

nn : =slze+1 ; Boxuser : luen : =f alse : one : = 1 ; 
CHARACTERISTIC EQUATIONCAI A , SIZE ,DEN ) : 

“ , (call polynom.inc) 
clrscr ; TextColor (yellow) ; , ,, , , 

writeln( ’ Denominator of Y(s)/U(s) :’)J 

vertpos : =2 ; 
for i : =nn downto 1 do 
begin 

3 : =nn-i ; 

g oscounter:=(j mod 4) +1; 

orizpos : =poscounter * 1 4 ; 
if poscounter = 1 then vertpos := vertpos+2; 
if i <> 1 then 
begin 

gotoXY(horizpos-1 0 , vertpos ) ; 
wrltetDENliJ :7:4) ; , , 

msg(’S +’ ,horizpos , vertpos ) ; str ( i-1 : 2 , exponent ) 

msg( exponent , hor izpos+1 ,vertpos-1 ); 
end 
else 
b 0gl n 

g otoXY ( horizpos-1 0 , vertpos ) ; write (DEN [i] : 7 : 4 ) ; 

; 

end; write In ;writeln; write In ;writeln ; 

root_f inder ( size , den , realpart 1 , imaginarypart 1 ,one); 
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writeln(’The poles of the Y(s)/U(s) are:’): 
writef’ REAL PART’); 

writelnf ’ IMAGINARY PART ’ ) ,* 

L:=1 to size do 



for 1 
begin 

writeC ’ _ _ . 

wr ite( realpart 1 [i] : 7 : 4 j ; _ wr ite( ’ 
wr i teln( imaginarypart 1 [ i J : 7 : 4 ) ; 
end ; 

keyread(key ,keyold,not_erased) ;clrscr ; 

for i:= 1 to, size do, 

PI 0[i .size] :=B[i , 1 ] ; 

do 

do 



+0 



) 



); 



for 0 0 *.= 2 to size 
begin 

for i:=1 to size 
begin 

for l;=i to size^do , ^ , 

P10[i,o]:= P10[i,j] +A1A[i,l] 



to size do 



end; 
end; 

for i : = 1 
begin , 

CC[i] :=0.0; 
for ^^:=1 to^s:^ze.c^o 



* P10[l,k]; 



P10Ld .i] 

to size do 



* C[1 ,J] + CC[i] ; 



CC 

end; 

for i : = 1 
begin 

m: =nn-i : 

if CC[m] <> 0.0 then goto 100; 
end; 

100: Text col or ( yellow ) ; 

wr itelnC ’ Numerator of Y(s)/U(s) :’); 

vertpos : =2 ; 

for i:= m downto 1 do 

begin 

J : = m - i ; 

poscounter : = ( d mod 4) +1 ;horizpos : =poscounter * 1 4 ; 
if poscounter = 1 then vertpos := vertpos+2; 
if i < > 1 then 
begin 

gotoXY ( hor izpos- 1 0,vertpos);write(CC[i]:7:4); 
msg(’S +’, hor izpos , vertpos ) ; str ( i-l : 2 , exponent ) ; 

msg( exponent ,horizpos+1 ,vertpos-1 ); 
end 
else 
begin 

g otoXY (horizpos-10,vertpos);write(cc[i]:7:4); 

; 

end; writeln; writ eln ;writeln;writeln; 
msize : =m- 1 ; 

root_f inder ( msize , CC , real root , imaginaryr oot , one ) ; 
wrltelnC’The zeros of the Y(s)/U(s) are:’); 



writeC ’ 
writelni ’ 
for i : = ) t 
begin 
write( ’ 



IMAGINARY PART ’ ) ; 
msize do 



REAl part ’ ) ; 



write(realroot[i] :7:4):_write( ’ 
wr iteln( imaginaryroot [i] : 7: 4 ) ; 



+d 



); 
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end; 

keyread(key ,keyold,not_erased) ;clrscr ; 
m : =m- 1 ; 

clr scr ; text col or (green ) ; highvideo ; 

msg( ’ Input the desired closed-loop characteristic 
, equation ’,1,6;; 

msg( ’ Factored (F) form or Coefficient (C) 
form ? ’ , 1 ,7) ; 

repeat read(kbd, input2 ) 

until inputs in [ ’ F ’ , ’ f ’ , ’ C ’ , ’ c ’ ] ; 

if (inputs = ’F’) or ( inputs* ’ f’ ) then 

begin 

if change then 

for j:=1 to size do 

begin 

str(real part [j ] : 1 0 : 2 , f ilvar [2^^ 1 +9] ) ; 
str ( imaginarypart [ j J : 1 0 : 2 , f ilvar [s»i+10] ); 
end; 

input_factored( ’POLES’ , size , r ealpart , 
goto 120 , lmaglnary_partT: 

end: 

if ( inputs* ’C’) or (inputs* ’o’) then 

clr scr ; m : =size+1 ; 

if change then 

for j : =m downto 1 do 

begin 

str(E[j] : 10:2 ,Filvar [size+22-j ] ) ; 
end ; umr an : = 1 ; 

input_coef f ( ’ POLES ’ , size , e ) ; 

root finder ( size ,E , r ealpart , imaginarypart , one ) ; 
clrs^r * 

write(^The roots of desired closed-loop 
characteristic ’ ) ; 
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writeln( ’ polynomial are: ’ ), 

write? REAL f>ART ’ ) ; 
writeln(’ IMAGINARY PART 

i:=i to size do 



for 
begin 
write( ’ 



;o size 



end; 



wr ite( realpart [i] :7:4):_write( ’ 
wr iteln( imSginary_part [i ] : 7 : 4 ) ; 



+j 



); 



); 



keyread(key ,keyold,not erased);clrscr; 
goto 180; “ 

120:for i:*1 to size do 
begin 

for,3:=1 to size do 

,j] :=0.0; 

end; 
end ; 

X ; = 0 • 

1 60 : i : *i+1 ; . _ 

if imaginary part [ij <>0.0 then goto 150; 



AA[i , i I : =real_part Li J ; 170; 

1 50 : AAfi , i J : = real part [ilj_ j:*i+1; 
-i , .1 J : =-imaginary_part Li J J 

[j , i] ; *imasinary_part [ij ; ^ 

*=.ll . AA[i , 1 J : =real part[ij; 
170:if (i-size) < 0 therl goto 160; 

CHARACTERISTIC EQUATION( AA , SIZE , E ) ; 



AA 
AA 
i 
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» H[j] ; 
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180:for i:=l to size do 
begin^ ^ ^ ^ ^ 

H[i]:= E[i] - DEN[i]; 
end; 

for i:=l to size do 
begin. _ _ _ 

5cc[i] :=C[1 ,i] ; 
end; 

INVERSE_FIND__MATRIX( PI 0 , CCC , SIZE , PIN , CCC , ER ) ; 
for_i:=l to size do 
:= 0 . 0 ; 

for 3* = 1 "to size do 
begin ^ ^ ^ ^ 

SHli] ;= HH[i] + PIN[j ,i] 
end; 
end; 

for i:=1 to size do 
begin 

for.3:=i to size do 

®|A?i,j]:= A1A[i,j] - B[i 
end; 
end; 

CHARACTERISTIC_EQUATION( AA , SIZE , E ) ; 

if CC[1]=0.0 then 
begin 

S ainl ;=1.0; goto 210; 
else 

galn1 : = E[1 ] / CCTl ] 

210:if gainl =0.0 the 
textcolor (yellow) ; 

to size do 



1] * HHCd]; 



hen gainl : = 1 . 0 ; 



for . i : = 1 

°®§H?i]:= HH[i] /gainl; 
end; 

140:for i:=l to size do 

begin Hj-iq /gainl . 

end;if HCsize];=0.0 then reduceorder : =reduceorder-l ; 
clrscr ;Writeln( ’Numerator of the Heq(s) is :’); 
vertpos : =2 ; 

for i:= size downto 1 do 
begin 

siz© “ij 

poscounter : = ( 3 mod 4) +1 ; hor izpos : =poscounter * 1 4 ; 
if poscounter = 1 then vertpos := vertpos+2; 
if i <> 1 then 
begin 

gotoXY(horizpos-l 0,vertpos);write(H[i]:7:4); 
msgC’S +’ .horizpos , vertpos ) ; str ( i-T : 2 , exponent ) ; 

msg( exponent , horizpos+1 , vertpos- 1 ) ; 
end 
else 
begin 

gotoXY ( hor izpos- 1 0 , vertpos );write(H[i] :7:4) ; 
end; 

end; writ eln ;writeln ; writ eln; 
reduceorder : =size-l ; 

root_f inder ( reduceorder ,H , realpart , 

, imaginarypart , one ) ; ^ 
writeln(’The roots of the Heq(s) are ;’); 
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write( ’ 
writeln( ’ 
for i : = 1 
begin 
writei ’ 



IMAGINARY PART ’ ) 
to reduceorder do 



Page 

REAL PART ’ ) ; 



write^ r^alpart [i2^7i^ ) f ’ 



Us 



writelnC imaginarypart [i] 
end; 

keyreadlkey ,keyold,not erased); clrscr; 
n2 : =NN-1 ; , “ 

190:if H[n2] <> 0.0 then goto 200; 
n2:=n2-1; goto 190; 

200 : n1 : =n2-1 ; 

TextColoriyellow) ; gotoXY(1,2); 

writelnC ’The feedback coefficients [ k J are 

writeln; for i;=1 to size do 

begin 

writeC ’ ’ ) ; write(HH[i] : 7 : 4 ) ; 

end;writeln; writeln; 



); 



); 



K' ] is : ’ ) ; 
) ; writeln ; 



writeC’The gain 
wr itelnCgainl ;7: 

keyreadCkey ,keyold,not erased); 
Repeat ~ 



5TDG3.X 

box printer rtextcolor (yellow) ;gotoxy(58, 15); 
writeC VP0L|.RES" ’); 
gotoxyC 1 ,25) ;write( ’ 
gotoxyU^,^); 

Read(Kbd, ch) ; 

If (chs’P’) or (ch=’p’) or (ch= 
then , Begin 

If (ch=’F’) or (ch=’f’) then 
Begin 

gotoxyC 45,15); t ext col or ( red) ; 

writeC ; PRINTING 

As sign (list, ’pole. res’ ); 



); 

F ’ ) or ( ch= ’ f ’ ) 



End 



{print pole nlacement output ) 

Rewr iteC list ) ; (to file on the current 



); 

drive ) 



If {ch=’P’) or (ch=’p’) then 
Begin 

gotoxyC 45,13); t ext col or ( red) ; 

writeC ’ PRINTING. .. , ’); 

Assign( list , ’ LST : ’ ) ; 

{print pole placement result to the printer) 

RewriteC list ) ; 

End; 

WritelnClist , ’ POLE _ PLACEMENT RESULT ’); 



xucxx.vxxoo. POLE PLACEMENT RESULT 

WritelnC list ) : WritelnC list ) ; 

WritelnC list .’ The Plant matrix A is : ’ 

writelnC list ) ; 

For i ;= 1 to size do 
Begin 

For 3 : = 1 to size do 

Begin 

Write(List , ’ ! , A1 A[i , 3 ] : 1 1 ) ; 

End; WritelnC list ) ; 

End; WritelnC list ) ; WritelnC list ) ; 

WritelnC list .’ The Input matrix 6 is :’); 
writelnC list ) ; 

For i := 1 to size do 
Begin 

For j ; = 1 to ni do 

Begin ^ ^ ^ 

WrlteCList . ; ’ B[i,d]:11); 

End;Writeln(list ) ; 

End; WritelnC list ) ; WritelnC list); 



); 
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Wr i teln( list. 
writelnC list ) 
For i 
begin 

for d : = 1 to 
begin 



= 1 to no 



Program Listing 
The Output Matrix 
do 
size do 

11 ); 



C is 



Write( list; j ! !»C[i,j] 



.St . ’ ’ 

end;wr iteln( list ) , 
end; WritelnC List );Writeln(list); 

wr iteC 1 ist Denominator of Y(s)/U(s) -’); 
writelnC list . ’ Descending powers of S:’); 
writelnC list ) ; 

=nn downto 1 do 



for 1 
begin 

wr iteC list . ’ 
end; writelnC lii 



,DEN[i] ^ 



Lvj., rrx j. a. j.st ); writelnC list y, , , , 

writelnC list The poles of the yCs)/UCs) 
writeClist,’ REAL PART’); 
writelnC list , ’ IMAGINARY PART ’); 

to size do 



are 



for i : = 1 
begin 
write! 
write! 
write! 
write! 



list, _ _ 

list , realpart 1 [i] ; 7 : 4 ) ; 

1 f ?t » ’ . . 

nC list , imaginarypai 



end; writelnC list ) ;writ elnC 



j;7:4;; 



7:4) 



wr iteC list Numerator of YCs)/UCs) -’); 
wr itelnC 1 ist . ’ Descending powers of S:’ 
writelnClist); 
for i : = m+1 downto 1 do 
begin 

wr iteC list , ’ 
end; writelnC list ) 



) 



writelnClist, ’The 
wr iteC list , ’ 



CC[i] :7:4) . , 
r itelnC 1 ist ) 



wr 
zeros 



of the 



YCs)/UCs) 

REAL 



are: ’ ) 
PART ’ ) 



do 



IMAGINARY PART ’ ) 



iteC list , ’ _ _ 

ite(list, real root [ij : 7 : 4 ) ; 

IteC list , ’ +3 ’ ) f T 

itelnC list, imaginary root [i] : 7 : 



) 



4) 



writelnClist, ’ 
for i:=1 to msize 
begin 

wr iteC list 
wr : 
wr : 

end: writeinT list ) ; writelnC I ist ) , 
writeC list ,* Desired closed-loop Characteristic 
^ polynomial ’); 

writelnClist,’- Descending powers of S :’); 

writelnClist); 

for i:= size+1 downto 1 do 

begin ^ , 

writeClist,’ ’,E[i]:7:4); ^ 

end; writelnC list); writelnClist); 

wr iteC list ,’ The roots of desired closed-loop 
^ characteristic’); 

writelnClist,’ polynomial are: ’): 

writeClist,’ r6aL PART’); 

writelnClist,’ IMAGINARY PART ’); 

for i:=1 to size do 
begin 



write 
wr i 
write 
write 



list 



te(list,real part [i] : 7 : 4 ) ; 
teC list , ’ +3 ’ ) L 

telnC list , imaginar.v part i 



) 



list, imaginary part i] : 
end; writelnC list); writelnClist); 



7:4) 
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Wr i te ( 1 
wr iteln 
wr iteln 
for i ; = 
begin 
write 
end; writ 
wr itelni 
write ■■ * 
writ 
for 
begin 
write 
write 
write 
write 
end;wr it 
wr iteln( 



ist Numerator of the Heq(s) 
(list,’ Descending powers of 
(list); 

size downto 1 do 



IS 

S 








). 



eln(list),r,x 

list, ’The roots of the Heq(s 
St , ’ REAL PAR 

list,’ IMAGINARY PART 

to reduceorder do 

list , ’ 

list , realpart [i ] :7 : 4 
list,’ +3 






); 



); 



■^3 ^ i * 

n( list , imaginarypart Lil ; 7 : 4 ) ; 
eln( list ) ; wr iteln ( list ) : 
list, ’The feedback coefficients 
^ are : ’ ) ; 

list); for 1:=1 to size do 



[ k ] 



wr i teln( 
begin 
wr ite( 
end; wr it 
wri te( li 

wr iteln(x.j..svxf^a.j.xxi . i * ~y j x \ 
writeln(list); write In (list); 



list,’ ’ ) ;write( list ,HH[i] : 7 : 4 ) ; 
eln( list ) ;writelnf list ) ; ^ 

St , ’ The gain f K ] is : ’); 

list ,galn1 : 7 : 4 ) : 



Close( list ) ; 

End; _ _ 

Until Ch in [’Q’,’q’ ]; 

Assign( SvsFile , ’ svs . COM’ ) ; 

(re-execute SVS main program) 
Execute( SvsFile ) ; 

End . 
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FILE: TYPEDEF.SYS 



i * * * * 
* 

* 

* * * * 
const 



typp 



var 



Program Listing 



Page 1 



TURBO GRAPHIX version 1.00A 
Type definition module * 

Copyright (C) 1985 by * 

BORLAND International * 

MaxWor ldsGlb=4 ; 

Maxorder = 9 1 
MaxWindowsGlb= 1 6; 

MaxPiesGlb = 1 0 : 

MaxPlotGlb=200; 

MaxBlocks =10; 

Str ingSizeGlb=80 ; 

HeaderSizeGlb= 1 0 : 

RamScreenGlb ; bool ean= true ; 

CharFile : str ing[Str ingSizeGlb] = ’ 4x6 . fon ’ ; 
MaxProcsGlb=27 ; 

MaxErrsGlb=7 ; _ _ 

Extension : Str ing[4] = ’ . svs ’ ; 

wrkstrings St ring [Str ingSizeGlb] ; 

WorldType= record 

x1 ,y1 ,x2,y2:real; 
end; 

WindowType= record 

x1 ,y1 ,x2 ,y2 : integer ; 
header rwrkstring; 
drawn , top : boolean ; 
size : integer ; 
end; 

worlds=array [1 . .MaxWor IdsGlbl of WorldType; 
windows=ar r ay [ 1 . .MaxWindowsGlb 1 of _ WindowType 
PlotAr r ay=ar ray . 1 . . MaxPlotGlb , i . . 2] of real; 
character =array .1..3] of byte 
Char Array = array [32.. 126] 

Pi eType= record 

area: real : 
text rwrkstring; 
end; _ 

PieAr ray=ar ray [1 . .MaxPiesGlb] 

" -Jo: 



26] of character 



BackgroundArray=arr ay [0 



,Y=array [(^.^7] 

* . 1 0, 1 . . 1 0l of 

. 1 0 . 1 . . 1 0] of 
.10] of real ; 



_ of PieType; 
of byte; 



of boolean 
real ; 
integer 



LineStyleArr 
Aryl s=Arr ay 
Ary2s=Ar r ay 
Ary3s=Array 
Str2=String 
Str4=Str ing 
Str5 = StringL ^ 

Str20 = Str ing .20. ; 

Str25 = string .25. ; 

Str40 = String .40. ; 

Str80 = String 80 
Str255=String[25$] ; 

XI WldGlb ,X2WldGlb , YIWldGlb ,Y2WldGlb , 

AxGlb , AyGlb , BxGlb , ByGlb : real ; 

X 1 Ref Gib , X2Ref Gib , Y1 Ref Gib , Y2Ref Gib : integer ; 
LinestyleGlb ,MaxWorldGlb ,MaxWindowGlb , 
WindowNdxGlb : integer : 

XI Gib ,X2Glb ,Y1 GlbTY2Glb : integer ; 

XTextGlb , YTextGlb , VStepGlb : integer ; 

PieGlb ,DirectModeulb , ClippingGlD , 

AxisGlb , HatchGlb rboolean; 

MessageGlb , BrkGlb , HeaderGlb , TopGlb , 
GrafModeGlb : boolean ; 

CntGlb , Co lor Gib : byte ; 

ErrCodeGlb :byte ; 
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Program Listing 



lyGlb : LineStyleAr ray ; 
ay [0 . . MaxPr ocsGlb J of '“WrkString; 
:ay [0 . .MaxErr sGlb J of "'WrkStr ing; 



LineStyleAr rayGll 
Error Proc : arr: 

ErrorCode : array, I 
PcGlb : string [40] ; 

AspectGlb : real ; 

Graf Base : integer ; 
worldrwor Ids ; 
windowl : windows; 

Char Set : CharArray ; 

D : Ary2s ; 

A1 A , B , C ; Aryl s ; 

Len, Space ,Drive : Str2 ; 

Size.ni , no, Block 1 ,Key ,Keyold: Integer 
Not erased, Finished, Exit , Inser ton , FI 



F1 0 : Boolean ; 



Ch : Char ; 

Template .Answer , Previous Answer : Str80 ; 

SvsFile :File ; _ _ ~ , , 

P ,Filvar : Arr ay [ 1 . . 35] of str40; {menu prompts) 
worldndxglb : integer ; 
escape , retriev, luen ; boolean; 
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FILE: BODEPLOT . INC Program Listing 
procedure Bodeplot ; 

label 1; {label dealer ation for goto statement) 



var 



Code , I , Count , NumberDecades , Star tDecade , 
EndDecade : integer ; 

Wf , Wo ,Wi ,DeltaW,Gain 



PlotArrayl , PlotArra.Y2 , MagPhase Array , 
FreqArray : PlotArray; 



ZMagn , Z Phase , PMagn , PPhase , Phase 

TempX,TempY 

temp 

OpenLoop 

.1 » d d .kkl ,m, 1 , cnpoles , sizezeros , one 
d€ 



Real ; 



real ; 
real ; 
char ; 
boolean ; 
integer ; 
ary4 ; 



dencoef f , numcoef f , cdencoef f 
realpartpole , imagpartpole , realpar tzero , 
imagpartzero , crealpartpole , cimagpartpole : ary3s ; 
PSI : aryls; 

function Log( X : real ): real ; ^ 

{computes the base-10 logarithm of X) 

Begin 

If X=0 then Log:=0 else 
Log := Ln(X)/Ln( 1 0) ; 

End; 

function Expon( Y , X : real ): real j ^ 

{computes Y raised to X power) 

Begin 

Expon := exp( X * (ln{Y))); 
end; 

begin 

ClrScr; boxuser ; one:=1: 

Characteristic equation( A1 A, size ,Dencoeff ) ; 
for i:=1 to si^e do 
begin. ^ ^ _ 

PSI [i,size]:=B[i,l]; 
end; 

f'oi* to size do 

begin 

for i:=1 to size do 
begin 

1 :=size-j j+1 ; 

Ick 1 * = i * 

PSI ti , j J : =Dencoef f [kkl ] * B[i,l]; 

for 1:=1 to size do 

besin 

PSiCi ,d] :=PSI[i , j] +A1A[i,l] ^ PSI[l,kk1]; 
end; 
end; 
end; 

for i:=1 to size do 

^§umcoeff [i] : =0 . 0 ; 
for j:=1 to size do 
begin 

Numcoeff [i] : =Numcoeff [i] + PSI[j,i] * C[1,j]; 
end; 
end; 

for i:=1 to size do 
begin 

nij — siz©“^1 “i * 

if numcoeff tm] <> 0.0 then goto 1; 
end; 

1 : sizezeros : =m- 1 : 

Clr scr ; TextColor ( lightblue ) ; 

wr iteln( ’ » * Bode Plotting Parameters 

TextColor (yellow ) ; 
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MsgC’Open (0) or Closed (C) Loop Plot?’, 5, 5); 

p 0^P 0 S."t 

input ( ’ A ’ , ’ ’ , 45 , 5,2, true , FI ,F1 0 ) ; 

{sets flag OpenLoop if) 
temp := copy ( answer , 1 , 1 ) ; 

^ ^ {user , selects the open) 

if not(temp in [’0’,’C’J) then beep{ 350 , 1 50 ) ; 

^ {loop option for plot > 

until temp in [’0’,’C’J: 

if (temp = ’0’ ) then OpenLoop := true 

, else OpenLoop := false; 

Msg( ’What is the first frequency to be 
, plotted?’,5,7K* 

Msg( ’( example : .01, 1, 100, etc . ) ’ , 1 0 , 8 ) ; 

Input ( ’N’ , ’ ’ ,50, 7,^, true, FI ,F10); 

Val ( answer , Wo , code ) ; {Wo is the first plotted freq) 



Msg( ’ Input number of decades do you want 
plotted? ’,5,10}; 

Input ( ’N’ , ’ ’ ,51 , 10, 2, true, FI ,F1 0) ; 

Val ( answer ,NumberDecades , code ) ; 

root_f inder ( sizezeros ,Numcoeff , realpartzero , 
imagpartzero , one ) ; 

root_f inder ( size , Dencoef f , realpartpole , 
imagpartpole , one ) ; 



? ain:= Numcoeff [sizezeros+1 ] ; 

^or.i:=1 to sizezeros+1 do 

^Sumcoeff [i] : = Numcoeff [i] /gain; 
end; 

boxuser ; 

for i:=1 to maxorder do CDenCoeff[i] := 0.0; 
for i:=1 to SizeZeros + 1 do ^ , 

CDenCoeff[iJ := Numcoeff [i] * gain; 

for i:=1 to Size,+ 1 do ^ ^ , 

CDenCoeffLi] CDenCoeff[i] + Dencoeff[iJ; 

if Size > SizeZeros then CNPOLES : =Size 
{NPoles should always be) 
else CNPOLES : =SizeZeros ; 

{greater, but to be safecompute new denominator roots) 

root finder ( Cnpoles , CDenCoef f .CReal Part Pole , 

~ CImagPart Pole , one ) ; 

StartDecade := trunc( Log( Wo ) ) ; , 

{compute linear scale to plot ) 
EndDecade := StartDecade + NumberDecades ; 

{log numbers. Also ^figure step) 

Wf : = Wo * Exponi 1 0 . 0 ,NumberDecades ) ; 

DeltaW ;= Expon( (Wf /Wo ) , 0 . 01 25 ) ; 

Wi := Wo; 

for Count := 1 to 81 do(do 81 iterations ... arbitrary ) 



Begin 



if OpenLoop 
begin 



then 



compute bode numbers if openloop 
and later if closed loop) 



) 
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ZMagn:=1.0; ZPhase:=0.0; PMagn:=1.0; PPhase:=0.0; 

for i := 1 to SizeZeros do 
(compute magn and phase of zeros for freq step) 
begin ^ , 

ZMagn : =ZMasn * Sqrt(Sqr (RealPartZero[l] )+ 

Sqr ( Wi-IinagPar tZer o [ I J T ) J 
if RealPartZero [i J = 0.0 then 
ZPhase : =ZPhase+pi/2 . 0 else; 
begin 

if r ealpartzero [i ] > 0.0 then 

ZPhase : =ZPhase - pi tarctanC (wi- 

imagpartzero [i J )/( -realpartzero [i J ) ) 

else 

ZPhase : =ZPhase 



end; 
end ; 



e:=ZPhase + arctan((wi- 
imagpartzero[i J )/ ( -realpartzero [i] ) ) ; 



■? 



be 



r i : = 1 to Size do 
compute magn and phase 

f in 
Ma^ 



if Rea 




of poles for freq step) 
salPartPoleCi] ) + 

Ihen 



PPhase:= PPhase+pi/2 . 0 else 
begin 

if realpartpole [i J > 0.0 then,, 

PPhase : =PPhase - pi +arctan( ( wi- 

imagpartpole [i] )/ ( -realpartpole [ij ) ) 

else 

PPhase := PPhase + arctan((wi- 

imagpartpole [i J )/( -realpartpole [i] ) ) ; 

end; 

end; 



PlotArrayl [Count , 1 ] := Log(Wi); 

(fill plotting matrix with magnitude values) 
PlotArrayl [Count , 2 J := 20*LogLGain* 

(ZMagn/PMagn) ) ; 



PlotAr ray2 [Count 
PlotArray2 [Count 




=Log(Wi);(fill p 
:= ( 1 80/pi )*(ZPn 



hase matri 
ase-PPhase 




(next stmt covers freq wrap-around) 
if PlotArray2 [Count , 2J > 0 then 

Plot Array 2 [count , 2 J : = Plot Array 2 [count , 2] -360 ; 



Wi : = Wi * DeltaW; (increment freq step) 

end 



, else 

(perform same steps as above if closed loop requested) 
begin 

ZMagn : = 1 . 0 ; Z Phase : = 0 . 0 ; PMagn : = 1 . 0 ; PPhase : = 0 . 0 ; 
for i := 1 to SizeZeros do 

begin , , ^ . 

ZMagn : =ZMas:n * Sqrt ( Sqr ( ReaiPartZero [ I] ) + 

Sqr(Wi-ImagPartZero[I]T ) i 

if RealPartZero[IJ =0.0 then 
ZPhase : =ZPhase+pi/2 . 0 else 
begin 

if Realpartzero[i] > 0.0 then 
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' r-i\\ 

•realpartzero [i] ) ) ; 



ZPhase : =ZPhase - pi +arctan( ( wi- 

imagpartzero [i J )/ ( -realpartzero [ij ) ) 

0 1 s s 

ZPhase : =ZPhase + arctan((wi- 
imagpartzero[iJ )/(-] 

end ; 
end ; 

for i := 1 to CNpoles do 
begin 

PMagn ; = PMagn * Sqrt( Sqr ( CReaiPartPole [I] )+ 

Sgr ( Wi-CIrna 2 PartPole[I J )} ; 
if CRealPartPoleClJ =0.0 then 
PPhase : =PPhase+pi/2 . 0 else 
begin 

if crealpartpole [i j > 0.0 then 

PPhase : =PPhase - pi+arctan((wi- 

cimagpartpole [i J )/ ( -crealpartpole [iJ ) ) 

else 

PPhase : =PPhase +arct^((wi- 

cimagpartpole [i J ) / ( -crealpartpole [ i] ) ) ; 



end; 
end ; 



PlotArra.yl [Count , 1 ] := Log(Wi); 

{fill plotting matrix, with magnitude values) 
PlotArrayllCount ,2] := 20*Log(Gain * 

r -1 / PMagn ) ) ; 

PlotArray2 [Count , 1 J ; =Log(Wi );{ fill phase matrix) 

PlotArray2 [Count , 2] := ( 1 80/pi )*( ZPhase-PPhase ) ; 



if 



{next stmt 
PlotAr ray2TCount , 2 
PlotArray2 [count , 2 



overs freq wrap-around) 

> 0 then , 

: = PlotAr ray 2 [count , 2] -360 ; 



Wi ;= Wi » DeltaW; 
end ; 
end; 

PlotBode( StartDecade .EndDecade , NumberDecades , 
PlotArray 1 , PlotAr ray2 , OpenLoop ) ; 
end; {bodeplot) 
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procedure PlotBode( StartDecade ,EndDecade , 

NumberDecades : Integer ; 

Plot Array 1 .Plot Array 2 : Plot Array ; 
OpenLoop : Boolean); 



const 

MagnArray: array[1..12J of 

char = (’M»,'^A’,’G’,’N’’I’,’T’,’U>,’D’, 
^ >E’ : ’ >d’,’B’5; 

PhasArray: array[1..12J of 

char = (’P’,^H’,’A‘,’S>’E’,’ ’,’d‘,’e’, 

> 2 ’’ y y y y >)• 

FreqArray: string[1^] = ’ ’ FRilQUEN6Y (rad/sec)’; 



var 1 , 0 ,n 
ch 

x1 ,x2 

Delta 

MagLabel 

PhsLabel 

DecLabel 

Titlel , 

Title2 

Dump Graph 

w 

quit 

list 



integer ; 
char ; 
integer ; 
realj _ 
string [ 
string C 
str ing[ 




string[80] ; 
Boolean ; 
real ; 
Boolean; 
text ; 



function Log( X ; real ): real ; 
Begin 
If X = ( 



End; 



!0 then Log:=0 else 
= Ln(X)/Ln(10); 



function Expon(Y,X ;real ) :real ; 

{computes Y raised to X power) 

Begin 

Expon := exp( X * (ln(Y))); 
end; 



Procedure PrintGraphData; 

{prints numbers to a file or printer) 

Begin 

LeaveGr aphic ; Clr scr ; 
repeat 

Textcolor (white ); gotoxy ( 20 , 1 0 ) ; 

writelnf’ PROGRAM OUTPUT OPTIONS ’); 

gotoxyC 20 , 1 3 ) ; 

writeln(’<P> Printer output ’); 

Textcolor (yell o 
writelnC ’ Ch 

Textcolor ( white 
writelnC’ <F> Li 
gotoxyC 20 ,16); 
writelnC’ on the current drive ’); 

gotoxyC 20 . 1 7 ) ; 

writelnr*<Q> Quit ^ ^ , ’); 

gotoxy (42 ,15); textcolor (yellow) ; wr ite( ’ "BODE . RES" ’ ) ; 



w;; ^ gotoxy C 2 W. 1 4; 
ecfc Your Printer! ^ 
): gotoxyC20 , 1 5 ) ; 

St to Pile name 



); 



); 




ch = 
begin^ 

if (ch = 
begin 



or (ch = 
then 



f ’ ) or ( ch 



) or ( ch = ’ f ’ ) then 



gotoxy ( 24,15); textcolor ( red) ; 
PRINTING, 



write( ’ 



’ P ’ ) or 



’ ); 
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assign(list , ’Bode. RES’ ) ; 
rewrites list ) ; 
end 
else 
begin 

gotoxv( 24,13); text col or ( red ) ; 
urite( ^ PRINTING 

assign( list , ’ LST ; ’ ) ; 
rewritet lis-f;) ; 
end; 



); 



Title1:=(’ w (r ad ) 

, Phase ( deg ) ’ ) ; 

Title2: = ( ’ — ---4i » ); 



Gain ( db ) 



wr iteln( list , ’ 

write In (list ) ;writeln( list ) ; 
writeln(list,Title1); 
writeln(list.Title2): 
writeln( list ) ; wr iteln( 1 ist ) ; 
• = 1 to §1 do 



BODE PLOT RESULT ’ ) ; 



1 ]); 



3); 



for i 
begin 

w := expon( 1 0 . 0 , PlotArrayl [i 
wr itelnC list , ’ ’^w:11:3, 

PlotArr ay 1 C i , 2] : 8 : 

PlotArray2 [i ,2] : 7 ; 
if i= 47 then 
begin 

write(list,chr( 12) ) ; 
write in (list); 
writeln(list,Title1 
wr itelni list . Title2 / , , 

wr iteln( list ) ; wr iteln( list ) ; 
end; 
end; 
end; 

until oh in C’Q’.’q’]; 

En*t © r* Gr* 3.ph. i c i 

Twhen finished printing, go back to graphics 
and display graph) 
swapscr een ; 
close( list ) ; 
end; 



)i 



mode 



Begin 

initgraphic ; 
Def ineWindow 
Def ineWindow 
Def ineWindow 
Def ineWor ld( 
Def ineWor ld( 
Def ineWor ld( 
SelectWor IdC 
SelectWindow 
Se t Backgr oun 
SelectWor IdC 
SelectWindow 
Dr awBor der ; 




SetLineStyle 
For I : = 1 to 
Dr awLine 



For J ; =0 to 



( 1 ) ; 

5 do {draw horizo 
{ StartDecade . -60+ 
EndDecade , -60+( 2 
NumberDecades-1 



do 



tal graph lines) 

20*n, 

’»!)); 



Begin 



{draw vertical logarithmic graph lines) 
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1 to 10 



Program Listing 
do 



For I:= 

Begin 

Delta: -StartDecade + (Log(I 
Drawline (Delta , -60, Delta, 60 
end; 

end; , ^ 

SelectWindowC 1 ) ; <y-cLXis titles) 

For I:= 1 to 12 do 

Begin , r -i ^ 

DrawText ( 5,55+6*1 , 1 ,MagnArray [ll )i 
DrawText ( 630 , 60+6*1 , 1 , PhasArray [ I J ) ; 

®nd; , , , 

DrawText (250 ,195,1 ,FreqArray ) ; (x-axis title > 
For I : = 0 to 6 do <y-axis scale label) 

®itr( 60-20*1 : 3 ,MagLabel ) : 

_ .... MagLabel); 



oTir ow — , riai 

DrawText ( 12,1 3+2l 



x.xcLrrj.eXt( ic., I , 

Str(0-60*I:4, PhsLabel ) ; 

DrawText ( 600 , 1 3+28*1 , i , PhsLabel ) ; 
end; 

For I:= 0 to NumberDecades do 

{label the logarithmic scale) 

Begin, 

Str ( Trunc (Start Decade )+I : 3 , DecLabel ) ; 
DrawText (36+C 570 div NumberDecades) *1, 
, 186, 1, DecLabel); 

DrawText ( 30+( 570 div NumberDecades) *1, 
1 90, 1 , ’ 1 0’ ) ; 
end: , ^ 

SetLineStyle( 0 ) ; 

~ ■ “): SelectWorldC 2 ) : 

otArray 1 ,1 , -81 ,0,1 , 0 ) ; 

, .plot the magnitude) 
SelectWorldC 3 ) : 

SelectWindowC 3 ) ; 



SelectWindowC 2 
DrawPolygonC P 



SetLineStylei , , , 

DrawPolygon( PlotArray2 ,1 , -81 ,0,1 , 0 ) ; 

(plot the, phase) 
copyscreen; (save screen to memory) 

repeat until keypressed; 
quit := false; 

repeat (call for graph options menu) 

if OpenLoop then 

Graph_MenuC ’ Open Loop Bode Plot ’, DumpGraph , quit ) 

0 1 S 0 

Graph MenuC ’Closed Loop Bode ’, DumpGraph , quit ) ; 

If DumpGFaph then PrintGraphData; , 

(dump numbers if desired) 
until quit: , , 

LeaveGraphic ; (leave graphics mode) 

end; 
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{* This program 
procedure box; 
var i : integer; 
begin 

Highvideo ; Text( 
wr ite(chr (218) 
for i := 1 to ' 



Program Listing 
draw line for menu box. 



> 



!o lor (yellow) ;gotoxy( 1,2); 

; (draw upper left cor 

'7 do (draw upper horizontal 1 



corner 
ine 



begin 

wr i . _ ^ . , 

end; write(chr 



write(ch;:( 1 ^ . 

writeC chr ( 192)) 




gotoxy(1,5); 
for i:=,2 to 78 
write( chr ( 1 96 ) 
write(chr(2l7 
for i:= 3 to ^ 
begin 

gotoxyC 1 , i ) ;write( chr ( 179)); 

gotoxyC 79 , i ) ;write( chr ( 179)); 
end; 

gotoxyC 1 , 23 ) ; 
for i ;= 1 to 78 do 
wr ite( chr (196)); 
writeln( 



( dr aw vertical lines) 



en 




a; 



input promt) 
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Procedure Fi_and_gamma( var T:real; 

var fitary Is; 
var gamma: ary Is 

label 1 ; 
var 

A2.A3,A4,A5 

i,3fk 

test value, step 



); 



: aryl s ; 

: integer ; 
: real ; 



begin 

test value :=1 .0e-08; 
step : = 1 . 0 ; 
for i:=1 to size do 
begin 

3 : = 1 to size do 



for 
begin 
FI .i , . 
FI .i ,1. 
A4 4 ,.1_ 
A2 .1 ,.i. 
A5[i ,3 ] 
end: 
end; 

1 : f or i : = 1 
begin 

for j : = 1 

begin 



= 0 . 0 ; 

= 1 . 0 : 

= AlAri, 



l1; 



= IT /2! 0 j ’ * FI [i ,3 ] ; 

= T FI 1 . 1,3 ] ; 



to 

to 



size 

size 



A3ri , .1 ] : = ( T/steD ) * A4[i,i]; 

FI 4,.-j. :=Flti,.1 J + A3IIJ1; 

-.-I. xoih, -i ^ ^ Jf ((step+1.0) * 

+^C T / (step+1.0)) * A3 [ 1 , 3 ]; 



A2 [i , ’3 ' t =A 2 [i 



( step+2 . 
A5[I ' 



rii*w 






end; 
end; 

for i:=1 to size do 
begin 

for 3:=1 to size do 
begin 

^4[i.3] 5=0.0; . 






: = 1 ^ to size do ^ ^ ^ 

i ,d] : =A4Ci 4 I + A1A[l,k] * A3[k,:]]; 



for 
A4 

end; 
end; 

step:=step+1 .0; 
for i:=1 to size do 
for i : =1 to ^size d( 
if abs(A3[i,3]) : 
then goto 1 ; 

for i:=1 to size do 
for .1:=1 to, ni do 
GAMMALI , j] :=0.0; 

for i:=1 to size do 
for j : = 1 to ni do 

for k:=l to, size do , , , , 

GAMMA[i,j] :=GAMMA[i,3] + A5[i,k] 



(test_value * abs (Fl[i,j])) 
{ Initialize ) 



end; 



(calculate GAMMA matrix ) 

* B[k,3] ; 
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Procedure Matr ix_reduction( Var dumy square : Ary 1 ; 

Var n , tn , rank : Integer ) ; 

Label 142;Label 145;Label 235:Label 170:Label 180; 
Label 162; Label 198;Label 200;Label 224;Label 226; 
Label 228; Label 234; 



VAR 



reduction matrix 
compar e_value , temp_value , sum 



start 



: Aryl ; 

: Integer; 

: Real ; 



Begin 

For i := 1 to n do 
For .1:=1tomdo 

r eduction_matr ix[i , 3 ] := dumy_square [i , j J ; 



compare value := 1.00E-8; 
Jjc :=1; 

Nmin : = n - 1 ; 



If nmin = 0 then goto 235; 

For i := 1 to nmin do 
Begin 

For 3 := 33 c to m do 

Begin 

Jj : =3 ; 

For k : = i to n do 
Begin 

K.i : =k ; ^ . 

If abs( reduction matrlxLk.JJ) >=compare value 
then GoTo 142;“ 

End; 

End ; 

GoTo 235; 

142:If kj = i then GoTo 145; 

For 3 I := il to m do 

Begin ^ ^ 

Temp value := reduction matr ix[.i , 3 1 j : _ _ 

reduUtion_matrix[l , j 1 ] : =reduction_matr ixLkj , j 1 J ; 

reduction matr lx[kJ , j 1 ] :=temp value; 

End; 

1 45 : temp_value := reductlon_matr ixLi . j d J 5 
For ik := 33 to m do _ . . r. t 

reductionmatr ix[l , 3 k] ;= reduction_matr ix[i , 3 k] 

/ temp_value ; 

If 1 01 then GoTo 180; 
s 1 3.r* 2 “ 2 * 

162: for l’:= start to n do 

Begin _ _ 

Temp value := reduction matr ixl 1 , 3 J j ; 

abs( reduction matr ix[l , 3 3 ' J ) <= compar e_value 



If , 

then GoTo 170; 

For 11 :=33 to m do _ 

reduction matr ix[l , 11 J : = 

1 Ufa * r* 



'temp_value 



reduction matrix[l,ll] 
reduction matr lx[i , 11] ; 



1 70 :End; 

GoTo 200; 

1 80 : imin : =1-1 ; 

For 11 :=1 to imin do 

Begin 

Temp value := reduction m 
If abs( reduction matrix 
then GoTo 198;~ 



^^rixCli , 33 ] ; 

[I 1 . 30 ]) <= compar e_value 
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Prograjn Listing 
to_M do 



For 111 ;= jn — 

reduction_matrixLli , Hi J : = 
reduction_matrix[li , Hi] -reduction_matr ix[i , IH] 



*temp value; 



= i + 1 



1 98 :End; start 
GoTo 162; 

200; If .•] .T = m then goto 235; 

^ Jdc ;= 3d +1 ; 

End; 

For jk ;= dd to m do 
Begin 

Temp value ;= reduction matrix[n, jk] ; 

lit abs( reduction matrix[n, Ik] ) <= compare value 

For .11 ;= dk to m do _ _ _ ^ 

reduction_matrixLn, d 1 ] ' =reduct ion_matr ix [n , d 1 J / 

temp_value ; 

GoTo 228; 

226 ; End; 

GoTo 23$; 

228;For lf;= 

Begin 

Temp value 

If Sbs( r 

then GoTo 234; 

For .11 ;=ddc to m 
sdu ■ ■ 



1 to nmin do 

;= reduction matr ix[lf . ] j c] ; 
Sbs( reduction_matr ix[lf , H c] ) <= 



compare_value 



re'duct idn_matr ix*^lf , d 1 ] 



=reduction matrix[lf ,d 1 ] 
-temp value * _ _ 

reduction matr ix[n , ] d c] ; 



234 ;End; 
235 ;Rank 
For i ; = 
Begin 
Sum : = 

For d 

Begin 



t = 0 * 

i to n do 

0 . 0 ; 

= 1 t o m do 



Tf abs( reduction matrix[i • 1 ] ) <= compare value 
then reduction matrix[i,.j] ;= 0.0; ~ 

dumy square[i,1] := reduction matr ix[i - .1 ] ; 

Sum ;= sum + abs (reduction matr ix[i , d ] ) ; 

End; 

If sum >= compare value then rank ;= rank +1 ; 
End; “ 

End; 

Procedure matrix multipHcatlon( VAR matr ix, matr ixl , 

matr ix2 ; Aryl s ; 

VAR L,M,N ; Integer); 

var mult matr ixl , result matrix, mult_matrix2 ; Aryls; 

k,d,i : Integer; 



to m do 
=1 to 1 do 



Begin 

For d : = 1 
begin 
For 1 
begin 

mult matrixl[i,d] : 
end; 

For k;=1 to n do 

mult matrix2[d,k] ; 

End; 

For 1 ;=1 to 1 do 



= matrix[l , d ] ; 

= matrixl Cd ,k] ; 
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=1 to n do 



= c 



Begin 
For .0 
Begin 

result matrix[i, 3 J 
For k:= 1 to m do 

result matrixri.jj : 
mult matr ixT [i , k] 
End; ~ 

End; 

For i : = 1 to 1 do 
Begin 

For 1 : = 1 _to n do 

matr ix2 [i , i ] := result 

End; 

End; 



. 0 ; 

result matrix[i,j] + 
mult_matrix2 [k , j ] ; 

matrix [i , d ] » 
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procedure Inverse find matrix(VAR AtAryls; 

~ var XD0T;ary4; 



label 
var 



var ORDER : integer ; 
var AINV: aryls; 
var X:ary4; 

var singular : integer ) ; 
5;label 6;label 12;label 16;label 20;label 51; 



B2 : aryls; 

comp .temper ary_value :real; 
n , i , 3 , k , m : integer 



begin 
n: =1 ; 
for i;=1 
begin 

ror j : = 1 
be — • “ 



to order do 
to order do 



flNYCi j] ; = .0,0; 
B2[i,3] :=A[i,3J ; 



end; 
end; 

for i : = 1 
begin 



to order 
to order 



do 



end; 

for i:=1 to order do 
begin 

comp ; =0 . 0 ; 

6*:if’(abs (B2[k,i]) 



-abs(comp)) <= 0.0 then goto 5; 



=B2[k,i] ; 



comp : 
n : =k ; 

* 

if ik-order) <= 0 then goto 6; 
if B2[n.iJ = 0.0 then goto 51; 
“ ■ '0 then goto 5l ; 

0 then goto 12; 
to order do 



if 
if 
for 
begin 

temporary value ;= B2[i,mJ; 

B2ri,mJ :^B21.n,m]; 

B2[n,m] :=temporar.^ ya 



lue ; 

value : = AlNV[i,mJ; 

J ; 



temporary value := AiNV 
AINV|.i,mJT= AINV[n,mJ; 
AINV[n,m]:= temper ary_ 
end; _ _ 

temporary value : =X [i J ; 

X [i J ; =X [n] ; 



value ; 



X n 



^ ^i=tempQrary value., 

l2:Xti] :=X[iT /~B2ti,i]^ 

al r- - n 



: =B2[i , ij ; 
sr do 

flNV[i -m] := _AINY[i,m] / temporary value; 
B2 [i , mj : =B2 [i ,mj / temper ary_valuF; 



temporary value 
for m:=1 to order do 
begin 



end; 
begin 



=1 to order do 




order 
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begin ^ 

SinV[ 3 ,n] := 

B2[o ,n] :=B2 
end; 

1 6 : end; 
end; 

goto 20; 

51 : singular : = 1 ; 

20 : end; 



AINVr-j ,n] - temporary value * 
AINYli nl * 

Cd ’ temporary_value * B2[i,n] 
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Procedure llnear_equation( var AA:ar.y2; 

var nw,m: integer ; 
var BBjX:ary3; 
var K: integer; 
var U : ary2 ) ; 

label 1;label 2;label 3;label 4:label 5;label 71; 
label 81;label l20;label 8;label 6 ;label 10;label 130; 
label 7; 



Var 

i s , i t , kk , n , k 1 ,11, mm , 1 , j 

a1 

id 

b1 , cl , ic , w 



Begin {linear equation) 
n : =nw ; 
mm; =m+1 ; 

for i:=1 to n do 
Begin ^ ^ ^ 

i1 [i , mm] : = BB [i] ; 



for 

A1 






to m do _ 

; = AA [ i , 3 ] ; 



: integer 
: ary2 ; 

; ary3 ; 

: real ; 



End; 
k j “ 1 j 

if (n-m) >= 0 then goto 1 ; 
it : = n+1 ; 
n: =m; 

for i;= it to m do 
Begin 

for :l:=1^to mm do 
Ai [i»d] 0-0; 



End; 
1 : f o 
ID 




=1 to m do 
= i ; 



2 : kk : =k+1 ; 
is : =k; 
it ; =k; ^ 

B1 :=abs(A1 [k,k]); 
for i : =k to n do 
Begin 

f or d ; =k to m do 

^®ff? abs(A1 [i,d] )-B1 ) < = 
is : =i ; 

bT ; =abs^1 [i ,d ] ) J 
3 :End; 

End: 

If(is-k) <= 0 then goto 4; 

for d •= k to mm do 

Begin 

cl 1= AI fis ,d ] ; ,. T 
AI ti§»d J •= AI [k,d] ; 

AI [k ,3 ] : = cl ; 

End; ^ 

4:li*( 
ic 

ID[it] :=lc 
for i : = 1 to n do 
Begin _ _ 

Cl :=A1 ri,it] ; ^ 

Al[i,ltJ := Al[i,k]; 



if(lt-k) <=0 then goto 5; 

CkJ^MiCit]; 



0 then goto 3; 
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A1 [i ,k] : =c1 ; 

5:ii A1[k,k] <> 0 then goto 71; 
kk : =k ; 
k : =k- 1 • 

for 1 : =kk to m do 
A1 [j ,dJ := -1 .0; 

goto 6: 

71 :if (k-n) >0 then goto 120; 
if_(k-n1 < 0 then goto 81; 

AILn.mm] :=A1[n,mm] / A1[n,nJ; 



. . . if (,k _ „ 

if_(k-n1 < 0 then go^o 8 
I [n ,mir * ' - ^ 

§ 0X0 7 ; 

1:for j : =kk to mm do 
Begin 

Hkfd] := Al[k,i] / Al[k,k]; 
or 1 : = kk to n do 



*abs(w)) >= 0 



for _ 

Begin ^ , 

w:=A1 [i,k] * A1 [k,j]; 

A1 [i ,1 J , : =A1 [i^.T J -w; 
if ( abs( A1 Ti , i] ) - 0.0001 
^then goto 8 ; 

„ Airi,dJ ;=0.0; 

8 : End; 

End: 

if (k-m) >0 then goto 120; 
if (k-m) = 0 then goto 6; 
k : =kk ; 
goto 2j 

6: for i:=kk to n do 

if A1[i,mm] <> 0 then goto 120; 

7:k1 :=k-1 ; 

for is:=1 to k1 do 

Begin 

i : =k-is ; 
ii : =i +1 : 

for it:=ii to k do 
Begin 

for .1 : =kk to mm do, , , ^ , 

Al[i,dJ :=A1[i,d] -Al[i,it] *Al[it,d]j 



End; 

for i : = 1 to m do 
Begin 

for d • = "I "to m do 
Begin 

if (ID[d]ri) <> 0 then goto 10; 
XfiJ :=A1[d,mml; 
if (k-m) = 0 then goto 10; 
for ^is : =kk, to m do 

U[i,is-k] := A1[d,is]; 



1 0 :End; 

End; 
k : =m-k ; 
goto 150; 

1 20 :Wr i teln ; 
WritelnC ’There 
130:End; (linear 



are no eauat ions ’); Delay( 2000 ) ; 
equation) 
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Procedure Polynomial_of_roots( var n: integer: 

var rr,ri:ary3s: 
var cf:ary4); 

{This program calculates the coefficients of a 
polynomial from its roots) 

label 1; label 3; label 111:label 222; label 2; 
label 5;label 21;label 50; 



var 



ick,nn,mp , i , ii . 1 ,m,mm,mmp 

sumreal ,prreal 1 , sumimag,prreal ,primag 



: integer ; 
: r eal ; 
;ary6; 



Begin 

nn;=n+1 ; 

CF [nn] : = 1 . 0 ; 

for m:=1 to n do 

Begin 

sumreal : =0 . 0 ; 
sumimag: =0.0; 

goto 2; _ _ 

1 : J[1 1 :=J[1] + 1 ; 

2:if (1-m) >0 then goto 50; 
if (1-m) = 0 then goto 5; 
mm : =m- 1 : 

for i;=l to mm do 
Begin 

ii : =- i+1 I - 
JLii] : ; 

End; 

5 : prr eal : = 1 . 0 ; 
pr real 1 : = 1 . 0 ; 
primag:= 0.0; 
for 1 : = 1 to m do 
Begin 

prreal ; =prreal 1 ; 
ick:=jrij ; 
eal T : = 



prr eal *RR [ick] - primax 
primag:= prr eal *RI [ickj +pr imag*I 



prr< 



*RI[ick] 
R[ickJ ; 



50; 
1 ; 



End; 

sumreal := sumr eal+prr eal 1 ; 
sumimag : = sumimag+pr imag ; 
for i;=1 to m do 
Be^in ^ 

if (JLl4“n+m-l) >0 then goto 
if (J[l]-n+m-l) <0 then goto 
End; 
mp : =n-m+l ; 
mmp ; =mp mod 2 ; 
if mmp <> 0 then goto 111; 

CF[mp]:= sumreal ;goto 222; 

11T: CF [mp] : =-sumreal ; 

222:End; 
goto 21 ; 

50 :Writeln; 

WritelnC ’ There is an error in polynomial . 

calculation from roots procedure’); 
Delay(2000) ; 

21 :mmp:= n mod 2; 
if mmp <> 0 then goto 3; 
for i:= n downto 1 do 
CF[i] : =-CF[iJ ; 



3 :End; 
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(Plot Nyquist is a routine to draw the Nyquist plot 
from the data generated in the Nyquist procedure. } 

Procedure Plot_Nyquist ( Star tDecade , EndDecade , 

NumberDecades : Integer ; 

Fr eqAr r ay , PlotAr ray 1 , 

Plotarray2 , 

MagPhaseArray : PlotAr ray : 

BigPic , OpenLoop : Boolean ) ; 



var 

Xint,Yint, i,.1,n,code ; 
Ymin,Ymax, Xmin,Xmax ; 
XLabel ,YLabel 
Titlel , Title2 ; 

DumpGraph , quit : 

w , XminL , YminL : 

Gr aphWi dthX , Gr aphWi dthY : 
Xexponent , Yexponent : 
XexpLabel , YexpLabel : 
GraphArray : 

List : 



[iJi; 



integer 
real ; 
String 
Stringy 
Boolean ; 
Real ; 

Real ; 
integer : 
string[3 ] ; 
Plot Array ; 
text : 



function Expon( Y , X ; real ); real j 

(computes Y raised to X power) 

Begin 

Expon := exp( X * (ln(Y))); 
end; 



Procedure Pr intGraphData; 

(dump data used 

Begin 

LeaveGraphic ; Clrscr ; 

repeat , ^ ^ ^ 

Textcolor (white ); gotoxy (_20 , 1 0 ) ; 
■ PROGRAM OUTPUT OPTIO 



to make graph to printer) 



# 



writeln 

writeln( ’ Check 

Textcolor ( white ) ; 
writeln(’<F> List to 
gotoxy ( 20 ,16); 
wr iteln( ’ on t 

gotoxyC 20 » 1 7) : 

writeln( ’ <Q> ^ Quit 



NS 

writeln(’<P> Printer output 
gotoxy(20 j 14) ; 

Your Printer ! ^ 

gotoxy(20 ,15); 

File name 

the current drive 



otoxyr42, 1 5) ; 

extcoior (yellow) ;write( ’ "NYQUIST . RES" ’ ) 
oto^yj28, 17 ); 



); 
’ ) 
’ ) 
’ ) 



got 
read 
If 



or ( ch = 
then 



^bd , ch 
ch = ’F 
ch = ’p 
begin 

if (ch = ’F’) or ( ch = ’f’) then 
begin 

gotoxy( 24,15); textcolor ( red) ; 

write( ; PRINTING 

assign(list, (Nyquist . RES ’ ); 
rewrite( list ) ; 
end 
else 
begin 

gotoxy ( 24,13); textcolor ( red) ; 
write( ^ PRIntInG 

as sign (list, ’ LST : ’ ) ; 
rewrite( list ) ; 



f ’ ) or ( ch = ’ P ’ ) or 



); 



); 
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end; 
Tit lei 



=( ’ w (rad) 
Phase (rad) 



Magnitude 

Xplot 



YPlot ’ ) 






Title2 
wr i teln _ 
wr iteln 
wr iteln 
wr iteln 
for i : = 
begin 
w : 
write 



> ) • 

RESULT ^ ) ; 



St,’ NYQUIST PLOT 

list ):writeln( list ); , 

list.Titlel); Writeln(list,Title2); 
list): wr iteln( list ) ; 

1 to 81 do 



47 



if i = 
begin 
wr i 
wr i 
wr i 
wr i 
wr i 
end; 
end; 
end; 

until ch in E’Q’,’q’]; 
EnterGr aphic ; 
swapscreen; 

close( list ) ; 
end; 



FreqArray[i ,1 ] ; 
ln( list j^w: 9:3,* _ , _ 

■ rray Li , 1 J : 1 1 ; 



.XST. , w: ■ . 

MagPhaseArray Li , 1 
MagPhaseArray , 
PlotArrayl Li , l . : 
PlotArrayl [i ,2J : 1 
then 



0 



teln( list 
te( list , chr (12)); 
telnT list , Titlel ); 
telni list ,Title2 ) ; ^ 

teln(list); wr iteln( list ) ; 



Begin 
if 

^®^P?3] 
P 'A'. 
P .5. 
P 6 



(Plot Nyquist) 
not(BigPic) then 



(prompt for window parameters) 






151 1N00503-01 01 01 
1 51 2N00504-01 01 01 
1 51 3N00505-01 01 01 
1 51 4N00506-01 01 01 



Clrscr; TextColor(LlghtBlue); 

writeln( ’ *** NYOUIST PLOTTING PARAMETERS ***>); 

textcolor(yellow) ; 

writelnC ’= = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =’); 



Text CO lor ( white ) ; 

>rdinates for 



the 



HlghVldeo ; 

wr iteln ;wr iteln; Textcolo 
msg( ’ Input Plotting Coor 
Nyquist Plot’71,6); 

TextColor (yellow) ; 
wr iteln; wr iteln; wr iteln; wr iteln; wr iteln; 
wr iteln( ’ X-Minimum: ’ 

’ X -Maximum : 



wr iteln 
wr iteln 
wr iteln 



* Y-Minlmum ; 
’ Y-Maximum : 



Input Handler ( ’N0306 
wr iteln :wr iteln; 
wrlteln( ’Any changes 
, (Y""/ N):^) 

Input ( ’A’ , ’ ’ ,45, 16, 
If (answer = 



’ , Escape ) ; 

to these parameters? 



Input Handler ( ■ 
^"^”"""^31 ,Xmin, 
4] , Xmax . 



Val( f il vaF 
Val ( f ilvar 



! 



true , FI , FI 0 ) 



then 
’ , Escape ) ; 
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Val( 

Val( 

end 

else 

begin 

Xmax 



Xmin 

Ymax 

Ymin 

end; 



f ilvar [5] , Ymin , code ) ; 
f ilvar [ 6 ] , Ymax , code ) ; 



: = 100 ; 

{set default 
: = - 100 ; 

: = 100 ; 

: = - 100 ; 



values for 



"big picture" plot} 



INITGRAPHIC; 

niceaxes(xmin , xmax, ymin, ymax, 
n t = 1 * 

for i;=l to 80 do 
if (ab 
( abs 
if no 



); 



if ( absCplotarrayl [i , 1 1 ) > Xmax) or 

( absi plotarrayl [1 ,2] ) > Ymax) then n:= n+1 ; 

■ ' then n:= n-1 ; 

{use 1 extra point beyond graph border) 

Dr awPolygon( PlotArray 1 ,n,-80,0,1 ,0); 

{draw graph on screen) 

n I = 1 * 

for i;=1 to 80 do _ 

if ( abs(plotarray2 [i , 1 ] ) > Xmax) or, 

( abs(ploxarray1 [ 1 > 2 J ) > Ymax) then 



n : - n+1 ; 

if n<>1 then n;= n-1; 

Setlinestyle (1); 

DrawPolygonC Ploiarray2 ,n , -80 ,0,1 , 9 ) ; 



Repeat until Keypressed; 

{Put option menu on screen) 

quit;= false; 
repeat 

if OpenLoop then 

Graph Menu(’Open Loop Nyquist Plot ’ , 
DumpGraph, quit ) 

0X30 

Graph Menu( ’Closed Loop Nyquist Plot’, 
DumpGraph , quit ) ; 

if DumpGraph then Pr intGraphData; 
until quit; 

LeaveGraphic ; 

end; 
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{Graph Menu provides a window on screen ajid offers the 
user options to make a title, print the graph, print 
the numbers from the graph, or quit and return to the 
menu . > 



Procedure Graph_Menu( Ti tleWindowName : STR25 ; 

var DumpGraphData, quit 



boolean ) ; 



var 

Line 1 



Line2 , Line3 
TitlePrompt ; 



str ing[40] 



Procedure 
begin 

Text Co lor (White ) ; 
Center(’*** Graph 



Title *** ’ ,1 ,2,80 ) ,* 



[S]i= 



1 01 0A04001 -01 01 00 ’ 
1 01 2A04002-01 01 00 ’ 
1 01 4A04003-01 01 00 ’ 



text CO lor (yellow) : 
msg( ’Line 1 : ’ , 1 , 1 0 
msg( ’Line 2: ’ , 1 , 12 
msg( ’ Line 3 : ’ ,1 ,14 

text CO lor (green ) ; 
msg(’Type your title 

Input_handler ( ’N0103 



for your graph . ’ , 6 , 20 ) 
, escape ) ; 



Linel : = 
Line2 : = 
Line3 : = 



copy* 

copy! 

copyi 



f ilvar 
f ilvar 
f ilvar 






40) 
40 ) 
40) 



end; 



Procedure ShowTitle; 

{makes title block and writes title to block) 
begin 

copyscreen: ^ ^ 

SetLineStyle( 0 ) ; 

DefineWindowf 3, 1 1 ,20,40,60); 
DefineWorld(3,0,0,40, 16 ); , , 

SelectWor ld( 3 ) ; SelectWlndow( 3 ) ; 

Def ineHeader { 3 , Tit leWindowName ) ; 



{puts header on box) 



SetBackgr ound( 0 ) ; 
SetHeaderOn;DrawBorder ; 

Dr awTextWf 1,4,1, Line 1 ) ; 
DrawTextWi 1 ,8,1 ,Line2): 
DrawTextW(_1 ,12,1 , Line3 ) ; 
SetBreakOrf; SetMessageOf f ; 

repeat 

read(kbd, oh) ; 
case ord{ch) of 
{allow user t 
72 ; MoveVer 
75 : MoveHor(-1 
77 : MoveHor 
80 : MoveVer 

end; 

until ord(ch)= 13 



p move t 


itl€ 


( -4 , true 


h ] 


( - 1 , true 




t 1 , true ) 


« J 

f \ 


( 4 , true ) 


• J 

9 \ 



box anywhere 
up arrow) 
left arrow) 
right arrow) 
down arrow) 



on screen) 



{freeze box and continue with <return> key) 



end; 

begin 
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DumpGraphDat 

selectscreen 



:= False; 
aetiineStyle( 0 ) ; 

- underlying screen and display menu box) 
Def ineWindow(4 , 11 ,20,55,90); 

Def ineWor ld( 4,0,0,20,20); 

SelectWor ld( 4 ) ; SelectWindow( 4 ) ; 

Def ineHeader ( 4 , ’Graph Options Menu 



copyscreen ; 

{ save 



) 



Set Header On ; ' SetBackgrouhd( 0 ) ; DrawBorder 
~ ■ ,4,1 , ’<P> Print 



DrawTextWC 1 



DrawTextW 
DrawTextW 
DrawTextW 
DrawTextW 
repeat 
Option ; 
case ch 
’ P’ 



’ T A ^ 1 
» 1 0 , ] 

1 



(display menu 



raph to the 
options ) 



printer ’ ) 



,13 

,17 



<T> Make Title to the Graph’); 
<N> Print Table of Numbers’): 
used to Generate Graph’); 
<Q> Continue to the Program’); 



{interpret user input) 



without 
.print to 



menu box) 
printer ) 



N 



of 

begin 

swapscreen ; 

(redisplays screen 
hardcopy ( false , 1 ) ; (i 
copyscreen ; 
ch : = ’ P ’ ; 

end; 

begin , , 

leavegraphic ; (leave graphics screen) 

TitlePrompt; (prompt for title) 

entergraphic ; (return to graphics mode) 

swapscreen: (bring back graph) 

ShowTitle ;( display title box on screen) 
copyscreen ;( save graph with title box) 
ch : = ’ T ’ ; 

end ; 
begin 

DumpGraphData := True; 

(sets boolean to cause numbers ) 

(to be printed) 
ch : = ’ N ’ ; 

end; 



end; 

until ch in 
if ch = ’Q’ 



end; 



[ ’ P ’ , X 

then ’ Qui-f; 



T 



else 



’N’ , ’Q’ ] ; 
:= true 
Quit := false; 
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procedure input_Factor ed( Zeros_or poles :str5; 

NFactors : integer ; var real_part , 
imaginary_par t : ary5s ) ; 

(Begin processing Factored form Input-Internal routine) 
var 

i,j : integer; 

test : real ; 

begin 



p 


[11- 




= ’ 0905N01 01 1 -0001 01 ’ 


p 


.12. 




= ’2505N01012-000101 ’ 


p 


.13. 




= ’ 0907N01 01 5-0001 01 ’ 


p 


.14. 




= ’ 2507N01 01 4-0001 01 ’ 


p 


.15. 




= ’ 0909N01 015-0001 01 ’ 


p 


1 6 




= ’ 2509N01 01 6-0001 01 ’ 


p 


-17. 




= ’091 1N01 017-0001 01 ’ 


p 


.18. 




= ’251 1N01018-000101 ’ 


p 


.19. 




= ’ 091 3N01 01 9-0001 01 ’ 


p 


.20. 




= ’ 251 3N01 020-000101 ’ 


p 


.21 . 




= ’ 091 5N01 021 -000101 ’ 


p 


.22. 




= ’ 251 5N01 022-0001 01 ’ 


p 


.23. 




= ’ 091 7N01 023-0001 01 ’ 


p 


.24. 




= ’ 251 7N01 024-0001 01 ’ 


p 


.25. 




= ’ 091 9N01 025-0001 01 ’ 


p 


.26. 




= ’ 251 9N01 026-0001 01 ’ 


p 


.27. 




= ’ 0921N01 027-0001 01 ’ 


p 


.28. 




= ’2521N01028-000101 ’ 


p 


.29. 




= ’ 0923N01 029-0001 01 ’ 


p 


.30. 




= ’ 2525N01 050-0001 01 ’ 



ClrScr; TextColor (White ) ; (write screen titles) 
gotoXYC 1 ,24) ; 

Invvideo( ’ Press <ESC> to change an entry’); 
Center(’*** Desired characteristic Polynomial 
Input***’ .1 ,1 .80); 

HighVideo ; 

wr iteln ; TextColor ( green ) ; 
if Zeros or Poles = ’ZEROS’ then 

wr iteln( ’ NUMERATOR Transfer Function 
Input — FACTORED Form’ ) 

else 

writeln( ’ Character istic equation 

Input — FACTORED Form’); 



HighVldeo:writeln;writeln; , ^ 

for ,1: = 1 to NFactors do (type prompt strings) 
begin 

wr iteln( ’ s = + j ’ ) ; 

wr iteln ; 
end; 



str ( ( NFactors *2+1 0 ) : 2 , strg ) ; 

specification := concat ( ’ N1 1 ’ , strg ) ; 

Input handler ( specif icat ion , escape); 

(call the input handler) 



for j:= 1 to NFactors do 

(compute the zero values from) 
begin (.input handler string Filvar) 

valCfilvar [2*5+9 T , test , code j ; 
if code = 0 then Real Par t [ j J : =test ; 

(val conversion successful if code) 
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val ( f ilvar [ 2 * j + 1 0 ] , test , code ) ; , 

<= 0 , error otherwise} 
if code = 0 then Imaginary_Part [ 3 J : =test ; 

end; 

end; {procedure Input_f actor ed) 

{Begin processing Coefficient form Input-Internal 
procedure ) 

Procedure Input Coeff( Zeros or Poles ;Str5; 

NCoeff : integer : 
var Coeff: ary4 ) ; 



var 

N 6 oeff_old 
test: real; 



, { counters ) 

integer; {holds old poly order if 
, changing order) 

{holds "val" results until validated) 



begin 

P[21 ] 

P[22: 



23. 

-24. 

- 25 - 

.26. 

.27. 

28. 

[29J 

.30. 



0406N01 021 -0001 01 ’ ; 

{ InputrHandler _descriptors for coeff) 



1806N01 022-000101 
3206N01 023-0001 01 
4606N01 024-0001 01 
0408N01 025-0001 01 
1808N01 026-0001 01 
3208N01 027-0001 01 
4608N01 028-0001 01 
041 0N01 029-0001 01 
181 0N01 030-000101 



{form input) 



NCoeff_old:= NCoeff; 

ClrScr: TextColor (White ) ; {print screen titles) 
gotoxy( 1,24); 

invvideoC ’ Press <ESC> to chajige an entry’); 
gotoxyC 1,1); 

writeln( ’ *** Desired characteristic 

Polynomial input ***’); 

TextColor ( Green) ; 

if Zeros or Poles = ’POLES’ then 



wri teln( ’ NUMERATOR Transfer Function 
Input — COEFFICIENT Form’) 

else 

^®^^write( ’COEFFICIENT Form ’); 

TextColor ( lightmagenta) ; 
if not luen then 

writeln( ’ Highest degree coefficient 
MUST be 1.0’); 



end; 

HighVideo ; 

if NCoeff > NCoeff_old then 
begin 

for .i:=1 to NCoeff old - NCoeff do 

for i:=NCoeff old -*• 1 downto 1 do 
Filvar [i-f-1 J := Filvar[iJ; 



end; 



if NCoeff < NCoeff old then 
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begin 

for i:= 1 to NCoeff old - NCoeff do 
for i:= 1 to^NCbieff old^+ l^do 
FilvarCi] := Filvar[i+1]; 

end; 



V©p*tPos t — 4* * 

for i:=NCoeff+l downto 1 do 
begin 

i:=NCoeff+l - i; 

PosCounter := (i mod 4) + 1; 
HorizPos := PosCounter * 14; 

If PosCounter = 1 then VertPos 



VertPos 4 - 2; 



if i<> 1 then{prompts for coeff input) 
begin 

msgC’s +’ jHorizPos .VertPos ) ; 
str( i-1 ; 2 , Exponent ) ; 

msg(Exponent ,HorizPos+l ,VertPos-l ); 
end; 

end; 



str ( ( 20+NCoeff+l ) : 2 , strg) ; (sets up and calls 

input handler) 

specification := concat ( ’ N21 ’ , strg) ; 

Input_handler ( specification , escape ) ; 

for i ; = NCoeff +1 downto 1 do 
begin 

val (Filvar [NCoef f +22-.1 1 , test , code ) ; 
if code = 0 then Coef f [J J ; =test ; 
end; 

end; 
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procedure character istic_equation( var A3:aryls; 

var N: integer; 
var C2 : ary4 ) ; 

{This procedures computes the characteristic equation 
polynomials coefficients by using the principle-minor 
method) 

label 10; label 20; label 30; label 4-0; label 50; 
label 60; label 15; label 70; label 80; label 90; 
label 100; label 1; label 2; label 5; 



var 



B3,A1 

integer vector 

D1 0 : array [1..280] of 

mm, ii,l,m,i,kk,nr,nc,k1 ,i1i, 
counter , ml m, even ,nn , k , 3 1 ,n 1 ,m 1 



: aryls; 
: ary6 ; 



W W O X j l±L I LXl f V yXXXX f xV f ^ 

temp_value , det , det correction, 
determinant old,determi 



real ; 

; integer ; 
determinant , value c2, value :real; 



begin 

nn : =n+1 ; 

for i;=1 to nn 

begin _ 

n2[i] ;=0.0; 



do 



do 



end; ^ 

C2[nn] : =1 .0; 
for m:=1 to n 
begin 

k : = 0 ; 1 : = 1 ; 

integer_vector [ 1 ] : = 1 ; 
goto 2; 

1: integer vector [IJ : =integer vector [lJ+1; 
2: if (1-m) > 0 then goto 90; 
if (1-m) = 0 then goto 5; 
mm: =m-1 ; 

for i : = 1 to mm do 
begin 

ii:=i+1; __ __ 

integer vector [ii J : =integer_vector [ 1 J+ 1 5 
end; 

5: for i;=1 to m do 
begin 

for kk:=1 to m do 
l^egin . 

nr:=integer vector 
nc:=integer vector j 
B3ti ,kk] :=A3[nr ,nc. 
end ; 
end; 
k: =k+1 ; 






counter := 0; 

for i 1 i ; = 1 to M do 

begin 

for 3 1 : = 1 to M do 

begin , _ 

Alli1i,3l]:= B3[i1i,3l]; 
end; 
end ; 

for i 1 i : = 1 to M do 
Begin 

k 1 : = i 1 i : 

30:if A1[k1,i1i] <> 0.0 then goto 10; 
k1 ; =k1 +1 • 

if*(k1-M) <= 0 then goto 30; 
goto 40; 
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10:if (ili-kl) > 0 then goto 40; 
if (ili-kl ) = 0 then goto 70; 
for m1m:=1 to M do 
Begin 

temp valuei =A1 [i 1 i ,m1 mj ; 

A1 [i 1 i ,m1 mj : = AT [k1 ,m1 mj ; 

A1 [k1 ,m1 mj ; = temp value; 

End; 

counter := counter+1 ; 

70:ii :=i1 i+1 ; 
if ii > M then goto 20; 
for m1m:=ii to M do 
Begin 

if A1[m1m.i1i] = 0.0 then goto 80; 
value : =A1 [ml m, i 1 ij A1[i1i,i1i]; 
for n1 : = i1 i io M do 
begin 

A1Lm1m,n1J:= A1[m1m,nl] - A1[i1i,n1] 



* value; 



to M do 
* A1[i1i,i1i]; 



end; 

80 : End ; 

20 : End; 
det : = 1 . 0 ; 
for i 1 i : = 1 
begin 

det : =det 
end; 

det correction: = exp( counter * 
detl5rminant old:=det correction 
even;= counter mod 2T 
if even <> 0 then goto 60; 
determinant := determinant old; 

g oto 50; 

0 : determinant : =-determinant old; 
goto 50; 

40 : determinant : =0 . 0 ; 



ln(1 .0) ) ; 
* det ; 



50 : D1 0 [k] : = determinant ; 
for i : = 1 to m do 
begin 

1 : =m-i+1 ; 

if (integer vector [1] -(n-m+1 ) ) 



> 0 then goto 90; 



if (integer vector [1] -(n-m+1 ) ) < 0 then goto 1; 
end; “ 



ml : =n-m+1 ; , , v x 

value _c2 : =exp( m * ln(1.0)); 
even:= m mod 2; 
for i : = 1 to k do 
begin 

if even = 0 then 
begin 

value c2 : =value_c2 ; 

goto 15; 

end; 

if even < > 0 then 
begin 

value_c2 ; =-1 . 0 * value_c2 ; 

end ; 

1 5 : C2 [ml ] : =C2 [ml ] + pi0[i'l * value_c2 
value_c2 : =exp( m * ln(1.0)); 
end; 



end;goto 100; 

90 :writeln( ’ Error in characteristic equation’); 
1 00 : end; 
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{ This procedure uses modified BARSTOW method to 
calculate the roots of the polynomial } 

PROCEDURE r oot_f inder ( var n: integer; 

var A:ary4; 
var u,v:ary3s; ^ 
var ir; integer); 



{ label decleration for the GOTO statement) 

label 3; label 4 ; label 7; label 9; label 10; label 13; 
label I5;label 19;label 20; label 23;label 30; 
label 32; label 33; label 34; label 36;iabel 49; 
label 50; label 52; label 53; label 76; label 72; label 73; 
label 75; label 76; label 81; label 82; label 100; 



var 

irev,i,nc,m,nl,np, 3 ,i1 : integer ; 
p , q , r , f , e , cbar , d, qp ,pp ; real ; 
H,B,C; array [1..21] of real; 



Begin 

irev;=ir; {talte given values) 
nc : =n+1 : 

for ii=l to nc^do 
Hti]:= A[i]; 



p:=0.0; {initialization) 
q : =0 . 0 ; 
r ; =0 . 0 I , 

3:if H[1] <> 0 then goto 4; 
nc : =nc-l ; 

UlncT ; =0 . 0 ; 

V [nc J : =0 . 0 ; 
for i:= 1 to nc^do 
H[i] ;= H[i + 1] ; 

goto 3; 

4:if (nc-1 ) = 0 then goto 100; 
if (nq-^) then goto 7; 




<> 0 then goto 9; 



r ; =-HM 
^oto 

f?if (a'6s(H[nc-1 ] /H[nc] )-abs(H[2] /H[1 ] ) ) > = 

then goto 19; 
ir ev ; =-irev ; 
m:=nc div 2; 
for i : = 1 to m do 
begin 

nl ; =nc+l -i ; 

Fi=Hrnl^ 



0.0 



end: 






:tnl] :=H 
■[i] :=F; 

CULL i 

if Q <> 0.0 then goto 13; 

p : =0 . 0 ; 

goto 15:, 

13:p:=D/q; 

q ; = 1 . 0/ q; 

15: if r = 0.0 then goto 19; 
r : = 1 . 0/r ; 
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1 ]U 



e-1 0 ; 
no. ; 

uHC, ; 

= 0 . 0 ; 
0 . 0 ; 



C 

BLnc+ 

C [nc+ ^ - 
np =nc-1 , 

20 for j : = 1 
begin 

for i 1 : = 
begin 



end 



to 1 000 do 
1 to np do 



i=nc-i ' 
riI:=H 


• 


+ R 


* Bl 


ii + 1 ] 


[ij :=B 


.i] 


+ R 


» C[ 


.i + 1 J 



dl UL « 

if Ubs(B[ll/H[11 )-e) <= 0.e 
if C[2J <> 0.0 then goto 23; 
r : =r +1 . 0 ; 

f oto 30; 

3:r :=r - B[1 ]/C[2] ; 

30: for i1 := 1 to np do 
begin 
i 



0 then goto 50; 



;=nc-i1 ; . 
ril:= Hti. 


] - P 


» B[ 


!i + 1 ] 


- q * 


Bl 


:i+2] 


[ij:= B[i. 


J - P 


* C[ 


.i + 1 ] 


- q * 


c| 


.i+2] 




tlj^n goto 32; 

> 0.0 then goto 34; 

<= 0.0 then goto 33; 

-e) > 0.0 then goto 34: 
-e ) <= 0.0 then goto 70; 

* C[4] ; 



q-= q * (q 

goto 49; , r 

36:p:= p + CBr2] 
q:= q + (-B[2J » 
49:end; 
e : =e * 10; 

5 Qto 20; 



en goto 36; 

+ 1 . 0 ); 



*^C[3] 
cbar + 



C[3 






then goto 52; 



50 : no : = nc-1 ; 

V rnc] : = 0.0; 
if irev >= 0 

U[nc] :=1 .0/r ; 
goto _53i 
52 : U [no] : =r ; 

53:for i:=l to^nc do 
HtiJ:= B[i+1]; 

0 : no : =nc-2 



if ifev 



qp: = 1 

PP:=P4iq 
goto 73; 
72:qp: = 

pp: =p/2 . 



>= 6 then goto 72; 

/q; 

q » 2.0); 



73:f:=pp * pp - qp ; 
if f >=0.0 th< 



- . =f «...en goto 75; 

U .nctl ] : = -pp ; 

U .ncj : =-pp ; 
V.nc±1]:=_sqrt(-f); 

V[ncJ : =-V[nc+1 J ; 
goto 76; 

75: if pp <> 0.0 then goto 81; 
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U[nc+Il:= -sqrt(f) 
goto^82; , 




nc+1 !] : =- 
nc+' 






-1 • =„qP?U[nc + 1 ] ; 
J : =0 . 0 j 

' ' Jf 



0 .0 , 

or , i : = 1 ^to 
[i]:= B[i+2 



c do 



* 



abs(pp ) 



goto 4; 
1 00 : end; 



sqrt(f ) ) ; 
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